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Chapter 1 
Introduction 


Currently we can consider as a rigorously established fact the existence of 
dynamical systems with a small number of degrees of freedom (A^ >2 — 
in autonomous case, > 1 — for non-autonomous systems) for which un¬ 
der certain conditions classical motion cannot be distinguished from random 
motion PEIEIE]. Typical features of these systems are nonlinearity and 
the absence of any external source of randomness. Thus, using such syn¬ 
onyms for the term ’’random” as ’’chaotic”, ’’stochastic”, ’’irregular”, we can 
state that there are nonlinear deterministic systems for which these notions 
express adequately internal fundamental properties that comprise important 
and interesting subjects for investigation. For the last 30-40 years examples 
of chaotic motion have been detected in every held of natural science, and 
their number continues to grow. 

The mechanism which is responsible for the existence of chaotic regimes 
in purely deterministic systems is based on local instability. It leads to ex¬ 
ponential divergence of initially close trajectories 

dit) = di0)e^^ 

where d is the distance between two points in phase space that belong to 
different trajectories. It could be shown [U] that local instability leads to 
mixing i.e. splitting of time correlations 

lim R{f, if) ~ 

t^OO 

for arbitrary functions in phase space f{z),ip{z) {z(t) = [qit),p(t)]). A prin¬ 
cipally new result that was established by N.Krylov consisted in the un¬ 
derstanding that the averaged increment of local instability (h) determines 
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the time of correlations splitting 


K = (h). 

This eqnation exactly establishes the connection between the dynamics of 
the system and its statistical properties. In other words, we mnst nnderstand 
stochasticity as a rise of statistical featnres in the system as a resnlt of local 
instability. 

Local instability of trajectories makes the problem of assignment of initial 
conditions in classical mechanics as fnndamental as the nncertainty princi¬ 
ple makes the problem of measnrement in qnantnm mechanics. From the 
traditional viewpoint hnite precision is limited by the imperfection of the 
measnring instrnment. Bnt the new viewpoint is different: we conld not 
record the resnlt of the measnrement with inhnite precision even if the ideal 
measuring instrument were constructed — there is not enough energy, time 
and paper. Uncertainty of measurement could not be completely eliminated 
even as a theoretical idea because inaccuracy in the data-line of a purely 
deterministic system is not produced by external randomness or noise but is 
connected with the hnite precision of initial conditions assignment instead. 
In particular, at the beginning of the last century Poincare showed [7] that 
for some astronomical systems a tiny imprecision in the initial conditions 
would grow in time at an enormous rate. Thus two nearly-indistinguishable 
sets of initial conditions for the same system would result in two vastly dif¬ 
ferent hnal predictions. In the presence of local instability this inaccuracy 
leads to very peculiar behavior of deterministic systems. For most physicists 
’’dynamical instability” and ’’chaos” became convertible terms. 

Fundamental progress in the understanding of classical nonlinear dynam¬ 
ics caused numerous attempts to integrate the concept of stochasticity into 
quantum mechanics. The core of the problem lies in the fact that energy spec¬ 
trum of every quantum system is discrete and thus motion is quasiperiodic 
if the system demonstrates hnite motion. On the other hand, the correspon¬ 
dence principle demands the possibility of transition to classical mechanics 
which demonstrates not only regular regimes but chaotic too. Several im¬ 
portant results that clarify this contradiction were established. However, 
long before complete solution of the problem there has been interest in its 
reduced variant ~ the search for peculiarities in the behavior of quantum sys¬ 
tems which have chaotic classical analogues. These peculiarities are called 
Quantum manifestations of classical stochasticity (QMCS). 
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Many decades would pass before dynamical chaos ideology was realized 
by the science community as a whole. According to the old ideology: 

• Chaos is an attribute of a compound system. 

• In any compound system it is possible to hnd out the elements of chaos. 

• Useful information is contained in those few places in which chaos is 
absent. 

• Physicists must look for non-chaos. 

The new ideology changed the situation principally: 

• Chaos is universal inalienable property of simple deterministic systems. 

• Chaotic dynamics is the most general way of evolution of an arbitrary 
nonlinear system 

• New interesting information is contained exactly in those branches of 
natural sciences where chaos is present. 

• Chaos is the major object of study. 

• The physicist must look for chaos! [8]. 

The basis of the present report is the new chaos ideology. In the context of 
this approach a general investigation of arbitrary nonlinear dynamical system 
involves the following steps: 

1. Investigation of the classical phase space, detection of chaotic regimes 
by numerical and analytical analysis of the classical equations of mo¬ 
tion. 

2. Analytical estimation of the critical energy for the onset of global 
stochasticity. 

3. Test for QMCS in the energy spectra, eigenfunctions and wave packet 
dynamics. 

4. Consideration of the interrelationship between stochastic dynamics and 
concrete physical effects. 
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The basic subject of the current report is to realize the outlined program 
for two-dimensional Hamiltonian systems with potential energy surface which 
has several local minima, i.e. multi-well potentials. 

Despite the huge number of papers concerning chaotic dynamics, Hamil¬ 
tonian systems with mnlti-well potentials have been somewhat neglected. 
Nevertheless the Hamiltonian system with multi-well potential energy snr- 
face (PES) represents a realistic model, describing the dynamics of transition 
between different equilibrinm states, inclnding snch important cases as chem¬ 
ical and nnclear reactions, nnclear hssion, phase transitions, string landscape 
and many others. Being substantially nonlinear, snch systems represent a 
potentially important object, both for the stndy of classic chaos and QMCS. 



Chapter 2 

Specifics of Classical Dynamics 
in Multi-Well Potentials — 
Mixed State 


Let us consider the characteristics of classical hnite motion in multi-well po¬ 
tentials. They are more complicated than in single-well potentials and allow 
the existence of several critical energies even for a hxed set of potential pa¬ 
rameters. This fact results in the so-called mixed state in such potentials [9]: 
at the same energy there are different dynamical regimes in different wells, 
either regular or chaotic. It is important to note that mixed state is a gen¬ 
eral feature of the Hamiltonians with nontrivial PES. For the hrst example, 
let us demonstrate the existence of mixed state for nuclear quadrupole os¬ 
cillations Hamiltonian [10]. It can be shown that, using only transformation 
properties of the interaction, the deformation potential of surface quadrupole 
oscillations of nuclei takes the form m- 

Uiao, as) = C^al + 2alra^iQal - air (2-1) 

m,n 

where Oq and as are internal coordinates of the nuclear surface during the 
quadrupole oscillations: 

= -Ro{l + ®oP2,o(^) y^) + as[L2,2(6', (fi) + Y2^-2{0, (^)]} 

Constants Cmn can be considered as phenomenological parameters. Restrict¬ 
ing to the terms of the fourth degree in the deformation and assuming the 
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equality of mass parameters for two independent directions, we get 
symmetric Hamiltonian: 


H = 


pI+pI 

2m 


+ UQo(x,y,a,b,c) 


( 2 . 2 ) 


where 


Uqo{x, y, a, b, c) = |(a;2 + + b (xy'^ “ y) + c{x^ + 

X = ao, y = \f2a2, a = 2 Cio, b = 3Coi, c = C20 

This potential is a generalization of the well-known Henon-Heiles potential 
[T^ . which became a traditional object for examination of new ideas and 
methods in investigations of stochasticity in Hamiltonian systems. It is es¬ 


sential that, in contrast to Henon-Heiles potential, motion in (2.2) is hnite 


for all energies, assuring the existence of stationary states in the quantum 
case. Hamiltonian (2.2) and corresponding equations of motion depend only 


on parameter W = b‘^/{ac), the unique dimensionless quantity we can build 
from parameters a, b, c. The same parameter determines the geometry of 
PES 


X 


2\2 


UQo{.x,y]W) =—Ax +y)+xy - — + [x +y) 


(2.3) 


2H^^“' ' " ^ 3 

Interval 0 < hP < 16 includes potentials with a single extremum — min¬ 
imum in the origin that corresponds to the spherically symmetric shape of 
the nucleus. In the interval hP > 16 PES of Uqo contains seven extrema: 
four minima (one central, placed in the origin, and three peripheral, which 
correspond to deformed states of nuclei) and three saddles, which separate 
the peripheral minima from the central one (Fig,2.1). Numerical calcula¬ 


tions of equations of motion in the region 0 < hP < 16 (region of single-well 
potentials) indicate regularity-chaos-regularity (R-C-R) transition: gradual 
transition from regular to chaotic motion when energy increases and restora¬ 


tion of regular motion for high energies (£g,2.2). In the next section we will 


discuss in details possible stochastization scenarios and methods of critical 
energy calculation. 

In the region hP > 16 (multi-well potentials) we confront with a substan¬ 
tially more complicated situation. In Figure |2]^ there are presented Poincare 
surfaces of section (PSS) for different energies. They demonstrate evolution 
of dynamics in central and peripheral minima of QO potential Uqo (with 
W = 18, when depths of central and peripheral minima are equal). At low 
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Figure 2.1: The level lines of the PES (2.3) for different structurally stable 
domains. 



Figure 2.2: The R-C-R transition in the quadrupole oscillations potential 
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energies motion is clearly quasi-periodic for both minima. Let us pay at¬ 
tention to the difference in topology of PSS. In the central minimum, PSS 
structure is complicated and has hxed points of different types, while in pe¬ 
ripheral minima PSS has trivial structure with only one elliptical hxed point. 

When energy increases, gradual transition to chaos is observed, but changes 
in character of motion are totally different in different minima. In the cen¬ 
tral minimum already at energy equal to a half of the saddle energy Es, a 
sizeable part of the trajectories is chaotic, and at the saddle energy there are 
almost no regular trajectories at all. At the same conditions, in the periph¬ 
eral minimum motion remains quasiperiodic. Furthermore, even at energies 
higher than the saddle energy there is a substantial part of the phase space 
occupied by quasiperiodic motion. In other words dynamics above the bar¬ 
rier has some kind of ’’memory”: the structure of phase space at energies 
greater than the saddle one is determined by the character of the motion in 
the local minima. 

Thus, in the energy region Es/2 < E < Es, classical dynamics is clearly 
chaotic in the central minimum and remains regular in peripheral minima. 
This type of dynamics, when chaoticity measured at hxed energy signihcantly 
differs in different local minima, represents the common case situation in 
multi-well potentials and is called the mixed state. 

As an example of the concrete realization of the (2.1) potential, let us 
consider deformational potential, that describes quadrupole oscillations of 
Krypton isotopes. Seiwert, Ramayya and Maruhn [13] restored the param¬ 
eters of the deformation potential of the quadrupole oscillations, including 
the sixth degree terms in deformation for isotopes /^7-74,76,78,80 


Uqo{x, y; a, b, c) = + y^) + b (^xy^ - t ) + 

-1-6 (^xy‘^ — yj + y‘^) + e (^xy‘^ — ^j -f f{x‘^ + y'^Y- 

The big experimental values of energy of the hrst 2’''-states for nuclei 
Ar74’76 indicate a spherical shape of nucleus surface, while the probabilities 
of the electromagnetic transitions 2~^ —)■ 2° and very low energies of hrst 
rotational states imply the possibility of superdeformation. The nonlinear 
ehects, which are connected with the geometry of PES must be exhibited in 
the superdeformed nuclei at relatively low energies of excitation. The PES 
of Krypton isotopes are presented in Fig j2.4[ 

As can be seen, inclusion of higher terms of expansion leads to consider¬ 
able complication of the PES geometry: for all considered Krypton isotopes 
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Figure 2.3: PSS for motion in the potential Uqq (2.3) with hF = 18 at 
different energies: E/Es = {0.5, 0.8,1, 2}. 
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Figure 2.4: The PES of Krypton isotopes. 
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Figure 2.5: Level lines and Poincare sections for and Dj potentials . 


we run into potentials with PES of complicated topology with many local 
minima. Its clear that, for all potential surfaces, C 3 .„-symmetry is preserved. 

The mixed state, which was shown above for the potential of quadrupole 
oscillations, is the representative state for a wide class of two-dimensional 
potentials with few local minima. According to the catastrophe theory na, 
a rather wide class of polynomial potentials with several local minima is 
covered by the germs of the lowest umbilical catastrophes, of type D^, 
D-j, subjected to certain perturbations. Let us note that the Henon-Heiles 
potential coincides with the elliptic umbilic . Fig j2.5| represents level lines 
and Poincare sections at different energies for multi-well potentials from a 
family of umbilical catastrophes and D-j-. 

Uds = 2|/^ — x'^ + xy^ + ^x'^ -I- 1 (2.4) 

Uor = -f- xy"^ - ^-x'^ + 

o Z v 

The mixed state is observed for all considered potentials of umbilical 
catastrophes in the interval of energies Ecr < E < Es (here Ecr is the critical 
energy of the transition to chaos in that local minimum, where chaos is 
observed at energies smaller than the saddle energy). 
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Chapter 3 

Chaoticity Criteria for 
Multi-Well Potential 

3.1 General formulation of the stochasticity 
criteria 

As noted above, stochasticity is understood as a rise of statistical properties 
in a purely deterministic system due to local instability. According to this 
dehnition values of parameters of a dynamical system, under which local in¬ 
stability arises, are identihed as regularity-chaos transition values. However, 
stochasticity criteria of such a type are not sufficient (their necessity poses a 
separate and complicated question), since loss of stability could lead to the 
transformation of one kind of regular motion to another. Despite this seri¬ 
ous limitation, stochastic criteria in combination with numerical experiments 
facilitate an analysis of motion and essentially extend efficacy of numerical 
calculations. 


3.2 Non-linear resonance overlap criterion 

First among the widely used stochasticity criteria is the nonlinear resonances 
overlap criterion presented by Chirikov [13]. The essence of this criterion is 
easier to explain by the example of a one-dimensional Hamiltonian system, 
which is subjected to periodic perturbation 

H = Ho(p, x) + Fx cos VLt. 
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For the unperturbed system we can always introduce the action-angle 
variables (/, 6) in which 

OO 

H = Ho{I)+ Xk{I)cos{ke -Qt), (3.1) 

k=—oo 


where 


27r 

xkii) = ^J deG^^x{i,e). 

0 


In the new variables the scenario of stochasticity, on which resonances 
overlap criterion is based, is the following. An external periodic in time 
held induces a dense set of resonances in the phase space of a nonlinear 
conservative Hamiltonian system. The positions of these resonances Ik are 
determined by the resonance condition between the eigenfreqnency 


u{I) 


dHo 

dl 


and the freqnency of the external pertnrbation For very weak the external 
helds the principle resonance zones remain isolated. As the amplitnde F of 
external held is raised, the widths Wk of the resonance zones increase 


Wk = 



i=ik 


and at F > resonances overlap. When this overlap occurs, i.e. under the 
condition 

+ hFfc+i) = \Ik — Ik+i\, (3-2) 

it is said that there is transition to global stochastic behavior in the cor¬ 
responding region of the phase space. The averaged motion of the system 
in the neighborhood of the nonlinear isolated resonance on the plane of the 
action-angle variables is similar to the particle behavior in the potential well. 
Isolated resonances correspond to isolated potential wells. The overlap of the 
resonances means, that the potential wells are close enongh to make possible 
the random walk of a particle between these wells. 

The ontlined scenario can easily be ’’corrected” for the description of 
the transition to chaos in the conservative system with several degrees of 



3.3. STOCHASTIC LAYER DESTRUCTION CRITERION 


19 


freedom. The condition of the resonance between the eigenfrequency and 
the frequency of the external held must be replaced by the condition of the 
resonance between the frequencies, which correspond to different degrees of 
freedom 


E 


dH„ 



0 . 


The intensity of the interaction between different degrees of freedom, i.e. the 
measure of nonlinearity of the original Hamiltonian, acts as the amplitude of 
the external held in this case. 


3.3 Stochastic layer destruction criterion 


This method could be modihed for systems with unique resonance [16]. In 
this case, the origin of the large-scale stochasticity is connected with the de¬ 
struction of the stochastic layer near the separatrix of the isolated resonance. 
The main point of modihcation consists of the approximate reduction of the 
initial Hamiltonian in the neighborhood of the resonance to the Hamiltonian 
of the periodically perturbed nonlinear oscillator 

X, t) = 2 ^^^ cosa: — Pcos k{x — r). 

The width of the resonance stochastic layer dehned by mi 



IT ~ 

M 

where 

2y/M 


P- nk ■ 

If P/M has an order then 

X 1 

IT ~ p e p 


where \ = 2k + 1 — s. And when 
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function W (p) has a bending flex point. Fast growth of W allows us to dehne 
the threshold of the stochastic layer destruction as 

Per = p + A — 1 j , 

when the tangent line in point pi to the function W[p] transects the p-axis. 

Application of these criteria in the presence of strong nonlinearity (which 
is inevitable when considering multi-well potentials) encounters an obstacle: 
action-angle variables effectively work only in the neighborhood of the local 
minimum. Because of this, interest in methods based on direct estimation 
of trajectories divergence speed, arises. One of the criteria of such a type 
is so-called negative curvature criterion [TH]. This criterion is based on the 
existence of the following scenario of the transition from regular to chaotic 
motion. 


3.4 Negative curvature criterion 

At low energies, the character of motion near the minimum of the potential 
energy, where the curvature is obviously positive, is periodic or quasiperiodic 
and is separated from the instability region by the zero curvature line. As the 
energy grows, the ’’particle” will stay for some time in the negative curvature 
region of the PES where initially close trajectories exponentially diverge. 
After a long time these result in a motion which imitates a random one and 
is usually called stochastic. According to this stochastization scenario, the 
critical energy of the transition to chaos E^r coincides with the lowest energy 
on the zero curvature line 

Ecr = Umin{K = 0 ) 


3.4.1 Results for one-well case 

Now lets demonstrate the efficiency of negative curvature criterion on the 
example of one-parametric family of potentials 


U(x, p\ p) = -(x^ + v^) + xy^ + 


(3.3) 


With p = —1 potential (3.3) reduces to Henon-Heiles potential and with 
p = -|-1 — to a separable one that is called anti-Henon-Heiles potential. In¬ 
teraction in any three-body system can be reduced to such a type of potential 
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Figure 3.1: Level lines and negative curvature region (in grey) for Henon- 
Heiles potential 


in the cubic approximation if its potential has the form 

U{xi,X2,X3) = U{xi- X2) TU{x 2- X3) + U(X3- Xi). 

Gaussian curvature of the considered PES turns to zero at the points that 
satisfy the condition 

+ (x - xof = Rl 

fl 


where 


_ 1 + h. p _ 1-F 

^0 — -^—) -^0 — 


4/i 


4/i 


When /i < 0, the zero-curvature line of (3.3) is an ellipse, which reduces to a 


circle for the Henon-Heiles potential (/r = —1) (Fig,3.1) 

2,2 ^ 

X +y = ^ 

On the zero-curvature line, energy is dehned by the expression 
U{K = 0) = ^fix^ + {1 + y,)x^ + ^{2 + y.)x+ 
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and the minimal value of energy on the zero-curvature line is 

= 0 ) = 

According to the considered stochastization scenario, this value offers the 
critical energy when the transition from regular to chaotic motion occurs. 
For the Henon-Heiles potential 

Umin{K = 0 ) = 

This result is in good agreement with numerical integration of the equations 
of motion, which is presented on Fig j3.2| in the form of the Poincare sections. 

It is necessary to make one important remark. Analysis of numerical 
integration of the equations of motion (e.g. in the case of PSS) allows us 
to introduce the critical energy of the transition to chaos, determined as 
the energy such that part of phase space with chaotic motion exceeds a 
certain arbitrary chosen value. Similar uncertainty is connected with the 
absence of the sharp transition to chaos when energy increases. Therefore 
a certain caution is required for comparison of the ’’approximate” critical 
energy obtained by numerical simulation, with the ’’exact” value obtained 
with the help of analytical estimations, i.e. on the base of different criteria 
of stochasticity. 

3.4.2 Failure for multi-well case 

Negative curvature criterion allows us to obtain a number of interesting re¬ 
sults in potentials with simple geometry (a single minimum) I2D1. However, 
when passing to multi-well potentials, this criterion fails to work correctly. 
In particular, for the above mentioned potentials and Dj), the structure 
of Gaussian curvature is similar in different wells. For example, for poten¬ 
tial according to negative curvature criterion we get the same value of critical 
energy for both minima Ecr ~ 5/9, but chaotic motion is only observed in 
the left well (see Fig j2.5D . 

A natural question immediately arises: is it possible, using only geomet¬ 
rical properties of PES but not numerically solving equations of motion, to 
formulate an algorithm for finding the critical energy for single local minima 
in multi-well potential? We will try to answer this question below in the 
framework of the so-called geometrical approach [231 EH- 
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3.5 Geometrical approaches 

The geometrical approach is based on application of differential geometry to 
study the chaotic dynamics of Hamiltonian systems. It turns out that as long 


as we consider the Hamiltonian of the form (2.2) we can restrict ourselves 
to the study of the conhguration space without loosing information. Thus 
application of this method for the analysis of the features of Hamiltonian 
dynamics in multi-well potentials seems natural because the characteristic of 
the multi-well potential is formulated in terms of conhgurational space. 

As is known, geodesics are among the main objects in Riemannian geom¬ 
etry. They are dehned as the shortest curves that connect two points on a 
manifold. The manifold in its turn is dehned by metrics. Having once hxed 
the metrics we thus dehne the distance on the manifold: 

= gikdx^dx^. 

We have the following condition for geodesics in this case: 


6 / ds = 0. 


(3.4) 


After variation we could obtain the differential equation for geodesics: 

= (3.5) 


ds ds 


where Tt^ are Christoffel symbols. 

Using variational principles it is possible to formulate Hamiltonian me¬ 
chanics in a geometrical way. Let us consider this question more closely. 
A trajectory of dynamical system is dehned according to Maupertuis 
principle: 


6 / 2Tdt = 0, 


(3.6) 


where 7 are all isoenergetic paths connecting end points, or according to 
Hamiltons principle: 
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To connect mechanics with Riemannian geometry we must choose the 
metrics that convert the expression under the integral into the length ele¬ 
ment. By such a procedure we will specify the manifold. Then trajectories 
will be geodesics on this manifold. We will call this CM - conhgurational 
manifold. This approach has an evident advantage: potential energy func¬ 
tion includes all information about the system, so one needs to consider only 
conhgurational space but not phase space. Lets emphasize that Christoffel 
symbols in this approach act as counterparts of forces in ordinary mechanics 
and metrics — as a potential. 

The simplest metric is the Jacobi metric. It has the form: 


gik = 2{E-Viq))5., 


ik 


(3.7) 


By means of this metric Maupertuis principle (3.6) could be rewritten in the 


form equivalent to the condition for geodesics (3.4) so that trajectories are 


dehned by equation (3.5). It could be shown that the geodesic equation (3.5) 


with Jacobi metric (3.7) leads to Newtons equations. 


Having equations of motion we now could consider local instability in 
geometrical form. Let q and q' be two nearby trajectories at f = 0. Then let 
us dehne the deviation 

Ji = (l'i- Qi- 

The dynamics of deviation are governed by the well known Jacobi-Levi-Civita 
(JLC) equation: 

dq^ kd(i 70 0^ 


ds"^ 


where is the curvature tensor. 


The two-dimensional case that we are 
interested in is very simple to consider because the only nonvanishing cur¬ 


vature tensor component is i?i 2 i 2 - To write JLC equation (3.8) explicitly we 


need to present local orthonormal basis. The simple choice is the following: 


ei = - 


dq^ dq^ 
ds ’ ds 


) ^2 — 


dq^ dq'^ 
ds ’ ds 


On this basis the deviation takes the form: 


J — ^ ^ ^i{s)ei(s) 
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It can be shown that in the two-dimensional case the JLC equation leads to 
two equations for deviation components on the chosen basis: 


+ 5^6 = 0 


^ — n 


(3.9) 


where R is scalar curvature. We can see that stability is determined by 
scalar curvature. For two-degrees-of-freedom systems Riemannian curvature 
has the form 

{E-V)AV + \AV\^ 

2{E-Vf 

where AV is positive for the considered potentials; thus Riemannian curva¬ 
ture is positive too. Due to this we cannot connect divergence of trajectories 
with negative Riemannian curvature. 

Pettini et ah 121 ] point out that instability is caused by oscillations 
of positive curvature and has parametric nature. Calculations of deviation 
dynamics of regular and chaotic trajectories in Henon-Heiles are presented 
in [22] • It is shown that, for regular trajectories, the normal component 
of deviation is bounded or grows linear in time, while for chaotic it grows 
exponentially. The initial conditions were chosen in specihc areas in the 
Poincare section. It is interesting to mention that initial conditions lying on 
the border of a regular island in section exhibit very slow exponential growth. 
This illustrates the so-called effect of sticky orbits. 


Let us briefly discuss the behavior of deviation in equation (3.9 second 


order for ^i). In order to do that we need hrstly to transform (3.9) to physical 
values, i.e. to time instead of interval izu. This procedure leads to the 
following equation: 


dR W dt ^ 


AI/ + 


(V^ 

w 


6 = 0 . 


(3.10) 


To simplify the equation and clarify its physical meaning, authors of paper 
make the substitution: 


(i(t) = Y(t)Vw ( 3 . 11 ) 

Obviously, 6(6 ^^0 Y(t) have the same behavior in the meaning of stability 


hence W{t) is bounded. Thus, using (3.11): 

(PY 


dE 


G{t)Y = 0 
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with 


n{t) = 


Ai/ + 


3(VK)^ 

2W 


4W^ 


dV dV 

qi -h q2 — 

qi q2 


1 2 


1 


E 

i,k 


d^V 

dqidq^. 


■qiqk 


Another geometrical approach, based on redehnition of covariant deriva¬ 
tive was proposed by Kocharyan in [22]. In the two-dimensional case these 
two approaches lead to the same equations for deviation. 

There are two possibilities for the solution to be unstable. The hrst 
case appears when Q{t) is negative. The second possibility is parametric 
resonance. As was mentioned before, the second case is more actual since 
the Riemannian curvature is positive. 

Fig|33| presents the distribution of {D{t)) ((...) means time averaged 
value) for an ensemble of 10^ trajectories in the chaotic well of potential 
with energy E = Es — 0.1. 

As we can see, there is no trajectory with negative {Q{t)). Further, only 
10% of trajectories have somewhere fl{t) < 0, but even for them {D,{t)) is 
positive. Although this is only a brief survey of the situation, we can expect 


that instability of the solution of (3.10) to have parametric nature. 


As was mentioned earlier, stochasticity criterion must derive critical en¬ 
ergy. However it was pointed out by Pettini et ah that we must know all in¬ 
formation about curvature oscillations to consider parametric instability. So, 
dynamical calculations are unavoidable in this investigation. Nevertheless, 
we could derive purely geometrical criterion by introduction of additional 
coordinates. Let us rewrite the JLC equation in the form which does not 
depend on dimensionality of the manifold: 


Id^llJl 

2 ds^ 




\jf- 


d 

ds 


= 0 


where is the so-called sectional curvature: 

E dq^ E dq"^ 

= —— 

and (J,v) = 0. Note that the point where iF^^^(J,v) < 0 is unstable. Since 
there are more than one sectional curvatures for the case N > 2, we could 
connect instability with the negative sign of some of them. It is assumed 
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Figure 3.3: The distribution of {kl{t)) for an ensemble of 10"^ trajectories in 
the chaotic well of potential for E = Eg — 0.1 










Regular and Chaotic Classical and Quantum Dynamics 


29 


that negativity of some of the sectional curvatures is sufficient condition for 
the rise of instability. One of the enlarged metrics is the Eisenhart metric: 


ds% = {gE)fiudq^dq'^ = Uijdq^dq^ — 2V{q){dq^)‘^ + 2dq^dq^^^ 


where a^- is the kinetic energy matrix. The additional coordinates are q^ = t 
and q^~^^ (the latter is connected with action). The nonvanishing components 
of curvature tensor are 


R 


OiOj 


d^V 
dq^dq^' 


Pettini et ah |23] considered on the constant energy surface 

for vectors 

J = (0, g2, -gi, 0), V = (1, gi, ga, gs)- 


Thus v) takes on the form: 




2iE-V{q)) \dql 


d^V 2 d^V 

ql + 


dq. 


2^1 - 


d'^V 

9gi9g; 


-gig2 


This value is easy to calculate at any point of phase space. In [ 23 ] the aver¬ 
aged value of v) is introduced. It is shown that for Henon-Heiles po¬ 

tential there exists a correlation between chaotic trajectories relative measure 
and averaged value of sectional curvature. This approach correctly predicts 
the value of critical energy. 

The case of multi-well potential is more complex. It is necessary to clar¬ 
ify whether this condition is sufficient for the development of chaoticity or 
not; clearly speaking we needs to answer the question, does the presence of 
negative curvature parts on CM always lead to chaos? Potentials with mixed 
state are a very convenient model for investigation of this question, since 
there exist both regimes of motion. So, we need to study, how the structure 
of v) differs in different wells. For that we calculate in [23| a part of 

phase space with negative curvature as a function of energy, i.e. a volume of 
phase space where v) < 0 referred to the total volume: 


Jdq<ipe(-A'W),i(g(q,p)-E) 

^ f dqdpd(H(q, p) — E) 

We carrie d ou t calculations for two potentials: and D-j. Calculations of 

ti{E) (Figj3.4[) show that there are parts, where < 0 in all wells. 
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Figure 3.4: Function fi{E) for (a) and Df (b) potentials. fi{E) for chaotic 
wells are represented by dotted lines, for regular — by triangles. 


but nevertheless chaos exists only in one well. Moreover, for the well with 
chaotic motion, function fi{E) gives the correct value of critical energy — at 
this energy fi{E) becomes positive. 

The situation with regular wells is more complicated. Although part of 
phase space, where iF^^^(J,v) < 0, is nonzero, chaos in the well does not 
exist. This can be viewed on the Poincare sections. For comparison in 
Fig |3.5| part of CS with negative Gaussian curvature (/ic) is shown. We can 
see that structure of negative Gaussian curvature is similar to the v)- 

structure. 

Investigation of the curvature of the manifold, as we can see from the 
above cited data, does not give a plain method for identihcation of chaos in 
any minimum, especially if there exist both regular and chaotic regimes of 
motion. It is impossible to determine a priori whether chaos existed in the 
system without using dynamical description (in our case which are Poincare 
sections). Nevertheless, we can efficiently use geometrical methods for in¬ 
vestigation of chaos in multi-well potentials. In thee potentials considered 
above chaos exists only in wells, which have two details: a non-zero part of 
negative curvature on the manifold and at least one hyperbolic point in the 
Poincare section. According to this, we can use the following method for 
identihcation of chaos and calculation of critical energy. At the hrst step the 
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Figure 3.5: Function ficiE) for (a) and Dj (b) potentials. 


Poincare section at low energy is drawn for the well and the presence of the 
hyperbolic point is detected. The quantity fi{E) must then be calculated 
(or the part of CS with negative Gaussian curvature). The value of energy 
at which n{E) becomes positive could be associated with critical energy. If 
there are no hyperbolic points in the section than chaos does not exist in the 
well. Consequently, geometrical methods could be efficiently used for deter¬ 
mination of critical energy in complex potentials and identihcation of chaos 
in general. However, we must carefully use these methods and combine them 
with qualitative methods, such as Poincare sectioning method. 


3.6 Normal forms 


The structure of the Poincare surfaces of section can be reproduced without 
resorting to the numerical solution of the equations of motion. For this, let us 
use the method of treating non separable classical systems that was originally 
developed by Birkhoff [2S] and later was extended by Gustavson [2B] . 

Every two-dimensional Hamiltonian near equilibrium point can be repre- 
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sented in polynomial form as follows 

^(p, q) = (p, q) + (q), 

H^^Kp,^) = EH{pl + ql), 

V{ci) = 

i<3 

q = p = {pi,P2)- 

The procedure of reducing to normal form depends on whether the fre¬ 
quencies Ui, are commensurable or not. If they are incommensurable then 
there exists a canonical transformation (q, p) —)■ rj) such that in variables 
{^,p) Hamiltonian r(,^,? 7 ) will be a function of only two combinations 

+ vl), ^ = 2 - 

In other words, the Birkhoff normal form is an expansion of the original 
Hamiltonian over two one-dimensional harmonic oscillators 

^(p, q) r(^, 7]) = uih + U 2 I 2 + E Oj^uI^ lu “1“ • • • 

If the frequencies u,y are commensurable, i.e. if there exist resonance relations 
of the type mui + nu 2 , the normal form becomes more complicated and will 
contain apart from other combinations of variables and rp, as well. Such 
extended normal forms are called the Birkhoff-Gustavson normal forms. 

We cite as an example the Birkhoff normal form (up to the terms of the 
sixth degree with respect to (C,h)) for the umbilical catastrophe in the 
neighborhood of right minimum: 


H{h, I 2 ) = 2 . 6 I 3/1 + 2 I 2 - 0 . 219 / 1/2 - O.OlT/f - 0.375/1- 
-0.005/f - 0.028/2/2 - 0 . 122 / 1/2 - 0.133/1 

The reduction of the Hamiltonian to the normal form solves the problem 
of the construction of a full set of approximate integrals of motion. The 
solution of the equations 


H{Px,Py,x,y) = E 

H{Px,Py,x,y) = /o 

X = const 
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E=1/2E^r 



a b 


Figure 3.6: The PSS for different values of energy for central minimum of 
obtained a) by numerical integration of equation of motion, b) by nor¬ 
mal forms. 


allows us to find the set of intersections of the phase trajectory with the 
selected plane {x = const) and to reconstruct the structure of PSS. 

The PSS for the quadrupole oscillations of nuclei which are con¬ 

structed in such a way, are shown in Fig j3.6[ The qualitative coincidence of 
topology of PSS calculated with the help of normal forms and the numerical 
integration of the equation of motion is noteworthy. 

Having failed in an attempt to formulate adequate criterion of stochastic- 
ity that is based on the estimation of the rate of divergence of the two initially 
close trajectories (nevertheless, we do not doubt the existence of such crite¬ 
rion), let us now return to the resonance overlap criterion. By means of it 
we will try to understand the differences in the phase space structure of the 
different local minima that realize different dynamical regimes. 
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As an example let us consider the application of this criterion to the 
dynamics that are generated by the potential of umbilical catastrophe D^. 
Accordingly the Hamiltonian in the reference frame connected with left (up¬ 
per sign) and right (lower sign) well has the form 


H = ^ix‘^ + + ^iy^ + ujjy^) + xy^ ^ V2x^ + 


1 

2 ' 


1 

4 ' 


where 


oji = ^2(2 =F V 2 ); U2 = 2. 

Now we make a canonical transform to the action-angle variables 


^ = V !7rCOS(pi; y=^^cosip2-, 
X = \/ 2 IiUi sin(pi; x = ^/ 2 TyJ^sm(p 2 . 


Thus Hamiltonian (3.12) takes the form 


H{Ii,l2,(pi,(p2) = Ho{Ii,l2) + fmim2 COs{miipi + m2(p2) 

ni,m2£A 

A: [0,1],[0,2],[0.3]. |0.4]. p.l],[2,-1] 


(3.12) 


(3.13) 


(3.14) 


where 


I 2 ) — Il(Xl + I2(X2 + ~ 1 



An item with indexes r = (ri, r 2 ) is called the resonance for a given value 
of energy E if there exist action variables (/[, J 2 ) such, that E = 
and 


nu^{ii,i;) + r2U2{ii,i;) = o 


where 


dHo 



1 , 2 . 


If the system is far enough from resonance, i.e. for all {mi, m 2 ) 


\miUi + miUi \ > /mima, 
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then avoiding the small denominator problem we could make a canonical 
transform to new action-angle variables thus eliminating angle dependence 
in lower orders under some small parameter. This procedure results in the 
redefinition of the integrable part Hq of the initial Hamiltonian and enlarge¬ 
ment of the set A of angle-dependent terms. After that we could encounter 
one of following three possibilities: 

1 . resonance terms are still absent in the considered region; 

2 . the single resonance term arises; 

3. multiple resonances arise. 

In the hrst case we should execute a new canonical transformation and keep 
carrying out the procedure until case 2 or 3 appears. In the second case 
critical energy of the transition to large-scale stochasticity could be dehned 
with stochastic layer destruction criterion [16]. And hnally in the third case 
we could use the Chirikov’s resonance overlap criterion [T5| to determine the 
critical energy. 

Small unbalancing near resonance 

\miUi +miUi\ < f^^rn2 


could be compensated by higher order terms that have arisen from the re¬ 
peated canonical transformation of the non-resonance terms. In the Hamil¬ 
tonian (3.14) the term [2, —1] is subjected to condition (3.13), and the con¬ 
sidered procedure leads to 


Hoih, h) — ~ 1 + i;:? ~ 


4a;i-|-5 


-py 


-1-3^2 _ 
2 


hh 

StJi 


In the left well let us take into account only two terms: resonance 

1/2 


1 / 2/2 , 

- COs2v7i-V72 

2oJi \0J2 


and ’’shaking” ones 


64a;i 


wi 1 


1 

CUi 


/lJ2COs(2(pi - 2(P2). 
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After this direct application of the stochastic layer destruction criterion leads 
to the value of the critical energy in the left well 


E ^ 

^cr — 


-0.51 


This value is in a good agreement with numerical integration results and in 
qualitative agreement with the estimation obtained by negative curvature 
criterion E^r = —5/9. Straightforward analysis of the integrable part of the 


Hamiltonian (3.14) shows that there are no resonance terms in the right 


well. Thus transition to large-scale stochasticity in the right well could be 
attained only when passing the saddle energy. This fact is in full agreement 
with numerical results. 

Resonance overlap criterion allows better understanding of the mechanism 
of the above mentioned regularity-chaos-regularity transition that exists in 
any (multi- or single-well) potential with a localized region of instability. Let 
us use similarity in structure of phase space of the considered two-dimensional 
autonomous Hamiltonian system with the compact region of negative Gaus¬ 


sian curvature and one-dimensional system with periodic perturbation (3.1) 

EH- 

The behavior of the widths of the resonances 


Wk = -{Wk+i + Wk), 


and the distances between them 


A4 = IT 


fc +1 


as a function of the resonance number is simplest when the satisfaction of res¬ 


onance overlap condition (3.2) for number ki (at a hxed level of the external 


perturbation) guarantees that this condition holds for arbitrary k > ki. This 
is precisely the situation, that prevails in the extensively studied systems of 
a Coulomb potential [ 30 ] and a square well [ 31 ] subjected in each case to 
a monochromatic perturbation. In the former case we have Wk ~ k^^^ and 
AJfc k~‘^/^, while in the latter we have Wk ~ k~^ and A4 [k{k + 1)]“^. 
As can be seen from Fig, 3.7 there is R-C transition (we will call this tran¬ 
sition a ’’normal” transition) for both the Coulomb problem and a square 
well, since there exists a unique point ki such that at k > ki the condi¬ 
tion Wk > A,Ik always holds. The motion is therefore chaotic. However, as 
the behavior of the widths of the resonances and of the distances between 
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Figure 3.7: The resonant spacing AR (solid lines) and the mean widths 
(dashed lines) as functions of the resonance numbers k. On the left: for one¬ 
dimensional Coulomb, on the right: for square-well potential. The critical 
point ki separates the regular range {R) from the chaotic one (S'). 


them as a function of the resonance number becomes more complicated, we 
can allow the appearance of an additional intersection point and thus a new 
transition: C-R transition, which we will call ” anomalous”. There is also the 
exotic possibility of the intermitting occurrence of the regular and chaotic 
regions in the phase space. 

We demonstrate that an anomalous C-R transition occurs in a simple 
Hamiltonian system: an anharmonic oscillator, subjected to a monochro¬ 
matic perturbation [321E7]. The dynamics of such a system is generated by 
the Hamiltonian 

H{p,xR) = Hq{j),x) + Fx cos Vtt (3.15) 

where the unperturbed Hamiltonian is 

p2 

Ho{p,x) = --h Hx" = E, {n = 21, I > 1) (3.16) 


The considered system fills a gap between two extremely important physical 
models: the harmonic oscillator (n = 2) and the square well {n = oo). 

In terms of action-angle variables the Hamiltonian Ho{p,x) (3.16) 
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becomes 123 


Ho{I) = 


2tt 


aG{n) 


G{n) = 


2yJ27rmT (1 + i) 

'4-r(i + a 


a = 


2n 

n + 2 


The resonant values of the action Ik that can be found from the conditions 


kuj{Ik) = uj{I) = 


dHo 

dl 


are 


Ik = a 


G{n) 

27r 


2n^ 


2na 

p 


, /5 = 


n — 2 


A classical analysis, based on the resonance overlap condition, leads to the 
following expression for the critical amplitude of the external perturbation 


per _ 2(2-3n)0 1^1 

^ 4:n /3 Xk 


G{n] 

271 


2n0 


^2n/9^4/3 ^^(n+2)/3 _ + l)(n+2)/3j ^ 


(3.17) 

where Xk is the Fourier component of the coordinate x{I,6). Expression 


(3.17) solves the problem of reconstructing the structure of the phase space 


for arbitrary values of the parameters. 

The phase diagram in Fig j3.8^ can be used to determine, at the hxed 
level of the external perturbation, the energy intervals of regular and chaotic 
motion. The snap-shot of E{x) at the right in Fig j3.8[3 conhrms that an 
anomalous C-R transition occurs. We can clearly see isolated nonlinear reso¬ 
nances which persist at large values of k, and near which the motion remains 
regular. The reason for this anomaly is explained by Fig.3.8. The plots 
of the resonance widths and the distances between resonances in this hgure 
demonstrate that there are two rather than one intersection points: k = ki 
and k = k 2 . Consequently, there is an anomalous C-R transition. 

Thus, for ID system with periodic perturbation R-C-R transition can be 
observed just as in the case of a 2D autonomous Hamiltonian system. The 
reason for the additional transition in both cases is common: the localized 
region of instability. In the hrst case this reason is the localized domain 
of overlap resonances, while in the second one this reason is the localized 
domain of negative Gaussian curvature. 
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Figure 3.8: a) The phase diagram of the R-C-R transition for Hamiltonian 
(3.15): resonance energy versus the critical values of the external pertur¬ 
bation F{n = 8). b) The snap-shot E{x) confirms that the anomalous C-R 
transition occurs. 
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Figure 3.9: The resonant spacings AI^ (solid line) and mean widths Wk 
(dashed line) as functions of the resonance numbers k{n = 8). There are two 
critical points ki,k 2 ■ The second critical point corresponds to the anomalous 
C-R transition. 





Chapter 4 

Quantum Manifestations of 
Classical Stochasticity — 
Formulation of the Problem 

After almost a hundred years of development, quantum mechanics became a 
universal picture of the world. On any observable scales of energy we could 
not hnd any violations of quantum mechanics. But this does not mean that 
from time to time quantum mechanics does not confront another challenge. 
The problem that arose at the face of quantum mechanics in the second part 
of the last century is called quantum chaos. The essence of the problem is 
the fact that, on the one hand, the energy spectrum of any quantum system 
with hnite motion is discrete and thus its evolution is quasi-periodic, but, on 
the other hand, the correspondence principle demands transition to classical 
mechanics which demonstrates not only regular solutions but chaotic too. 
This deep and serious problem requires an answer hrstly to the question: 
what does it mean when one theory is a limiting case of another? [28] 

Usually a more general theory G is connected with a special theory S 
with a dimensionless parameter 5 such as 

G 6^ Q 

For example, if under G we understand special relativity theory, and under 
S classical mechanics, then b jc?. In the simplest case we can represent 
the general theory as a Tailor series over parameter 5. 

However such a simple situation is a very rare exception. In the most 
general (and the most interesting) case the limit limG is singular and the 

5—^0 
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transition G —)• S' is far from being trivial. So, for example, the transition 
from Navier-Stocks equations (viscous fluid) to Euler equations (ideal fluid) 
is singular: dissipation does not turn to zero smoothly at zero viscosity. So 
difficult for investigation, the no man’s land between the two theories contains 
new physics, like turbulence or critical behavior at phase transitions. It is a 
similar region where we have to study the influence of classical stochasticity 
on quasi-classical behavior. In our case G stands for quantum mechanics, S 
— classical mechanics and 6 is some dimensionless combination of physical 
quantities with h in the numerator. According to Berry [2H] , the h —)■ 0 limit 
swarms with non-analyticities. 

As we are interested in absolutely concrete aspects of the semiclassical 
limit, namely: how the presence of classical chaos is reflected in quantum 
quantities, we shall discuss one more principal difficulty. As in classical me¬ 
chanics chaos is realized only on large time scales (required for complete 
mixing, i.e. for realization of the limiting tendency to zero of the correlation 
function), any useful discussion of semiclassical limit must simultaneously 
account for both the limit t —)■ oo and the limit h —)■ 0 |2H|. A natural ques¬ 
tion arises whether the two individually non-trivial limits t ^ oo and h —)■ 0 
commute? The answer is negative: long-time semiclassical evolution funda¬ 
mentally differs from long-time classical evolution — in the common case 
situation the classical long-time limit is chaotic, while in semiclassics the 
temporal asymptote is not, and any chaos represents just a transition pro¬ 
cess. Therefore in an attempt to construct the quantum theory of dynamical 
chaos we immediately confront a row of evident and very deep contradictions 
between well-established principles of classical chaos and fundamental prin¬ 
ciples of quantum mechanics. What is the reason for those contradictions? 

As is well known, the energy spectrum of any quantum system that un¬ 
dergoes hnite motion, is always discrete. And it is not a property of a con¬ 
crete equation, but a consequence of the fundamental principles of quantum 
mechanics: the discrete nature of the phase space or, more formally, the 
non-commutativity of quantum phase space. Indeed, according to the un¬ 
certainty principle, an individual quantum state cannot occupy the phase 
volume Vi < , where N is the dimensionality of the conhguration space. 

Therefore a motion limited by a region V will contain VjVx eigenstates. Ac¬ 
cording to existing ergodic theory such motion is considered as regular, in 
contrast to chaotic motion with continuous spectrum and exponential insta¬ 
bility. The latter statement can be verihed using the notion of algorithmic 
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complexity, which can be defined as the relation: 


C = 




N, 


(4.1) 


out 


where and Nout are expressed in bits input and output routine lengths 
respectively. This quantity can be determined for any moment of time; how¬ 
ever the distinction between regular and chaotic motion manifests only in the 
limit t —)■ oo. If the motion is chaotic, then C —)■ const > 0, if it is regular — 
C* —)■ 0. In order to understand the reason for that we should note that the 
output length of routine — information about the orbit in an arbitrary mo¬ 
ment of time — grows proportional to t. The input data sequence consists of 
two main parts. The first is the algorithm for solution of equations of motion, 
its length does not depend on t. The second is the definition of initial condi¬ 
tions with the precision required to reproduce the required final result. For 
chaotic systems, where errors grow exponentially, this part is proportional 
to t and therefore dominates in the input routine length. Therefore in that 
case the algorithmic complexity (4.1) will tend to constant. For non-chaotic 


systems part of the input routine, connected with initial conditions, grows 
slower (for example as Inf when the errors grow linearly), and the limiting 
value of the algorithmic complexity is C = 0. 

All experiments performed up to the present time showed strict fulfillment 
of that rule for classical chaotic systems. For quantum systems, as for those 
that are chaotic in the classical limit, and for those that are regular, only zero 
algorithmic complexity was observed. This result can be briefly formulated 
in the spirit of Bohr complementarity: classical evolution is deterministic, 
but random, quantum evolution is not deterministic and it is not random. 
In other words, the problem consists in the fact that the discrete nature of 
the spectrum never implies chaos, or more exactly any resemblance to chaos 
in the sense of the ergodicity theory, in any quantum system with finite 
motion. Meanwhile the correspondence principle in the semiclassical limit 
requires the presence of chaos, connected with the nature of motion in the 
classical case. 

If to state the point of view that chaos never appears in quantum me¬ 
chanics, then a possible reaction to that is just to give up the study of the 
question. But it will mean that we avoid the challenge that Nature gives 
us in the limit of small h and large t, which is equivalent to ignoring other 
singular phenomena, such as turbulence or phase transitions. An alternative 
point of view consists in the fact that not waiting for the complete solution 
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of the problem (or rather for its correct formulation) we can study its lim¬ 
ited variant: investigation of special features of quantum systems behavior 
whose classical analogues are chaotic , or, in other words, search for quantum 
manifestations of classical stochasticity (QMCS). It is the problem that will 
be considered in the following chapters on quantum systems with potential 
energy surfaces of non-trivial topology. 

Deterministic chaos is a general feature of Hamiltonian systems with the 
number of simple integrals of motion less than the number of degrees of 
freedom. Lack of the full set of integrals of motion (full set consists of such 
number of integrals that is equal to the number of degrees of freedom of 
the quantized system) makes the traditional procedure of multi-dimensional 
systems quantization unrealizable. Let us consider this statement in details. 

As is well known [33] , in the one-dimensional case it is always possible to 
introduce such canonically conjugate ’’action-angle” variables that the Hamil¬ 
tonian becomes a function of action variable only. The standard dehnition 
of the action variable concerns integral along the periodic orbit 



where p{x) is the particle’s momentum. In the context of the semi-cIassical 
approach we could construct an approximate solution of Schrodinger equation 
in the terms of the integral along classical trajectory [31| 



= —e ^0 

Vp 


(4.2) 


This construction makes sense only in the case when phase grows on multiples 
of 271 along the periodic orbit. This limitation immediately leads to the semi- 
classical quantization condition 



where n is nonnegative integer number and /i is the so-called Maslov index 
which is equal to the number of points along periodic orbit where the semi- 
classical approximation is violated (in the one-dimensional case this occurs in 
turning points and fi = 2). Semi-classical energy eigenvalues En are obtained 
by the computation of the Hamiltonian H{I) for quantized actions 
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For multi-dimensional systems this procedure could be executed only in the 
case when the number of integrals of motion is equal to number of degrees 
of freedom, i.e. for integrable systems. In this case the procedure is called 
the Einstein-Brillouin-Keller quantization. Let us restrict with the two- 
dimensional case for simplicity. If the system is integrable, then there exist 
two couples of canonically conjugate action-angle variables (Ji, ^i) and {I 2 , O 2 ) 
such that classical the Hamiltonian depends only on action variables 


// = //(/„/2 


(4.3) 


Finite classical motion is periodic in every angle variable with frequencies 

In the general case frequencies are not close to each other, and motion 
in four-dimensional phase space is quasi-periodic. Phase trajectories he on 
invariant tori that are dehned by integrals of motion R. Semi-classical wave 


functions could be constructed in the form that is analogous to (4.2), but 


turning points must be replaced by caustic surfaces. Uniqueness of the wave 
function demands the quantization condition 


Ii= in 




(4.4) 


As in the one-dimensional case energy eigenvalues could be obtained by 


the substitution of (4.4) into Hamiltonian (4.3). It was understood by Ein¬ 


stein in 1917 that this method could be applied only to quantum integrable 
systems with trajectories lying on tori. For non-integrable (i.e. chaotic) sys¬ 
tems consistent quantization method did not exist for half of the century. 
But how to implement quantization in the non-integrable situation? 

Progress in the problem of chaotic systems quantization was obtained 
with the help of Feynman’s formulation of quantum mechanics. The hrst 
mention of the applicability of path integrals to chaotic systems was given 
by Selberg (1956), who built the dynamics of the particle on the Rieman- 
nian surface with negative curvature in the terms of path integrals. This is 
certainly a chaotic system although this term did not exist at that time. 

Gutzwiller was the hrst who successfully applied an analogous approach 
to quantization of chaotic systems. In 1982, he had shown that semi-classical 
approximation in the form of path integrals allows us to obtain the spectrum 
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of a chaotic system |1]. This study was the culmination of the large series 
of his works [35l|36l|37lEHll39|. Works of Balian, Block [401 ITT] concern the 
same approach — they connect classical periodic orbits with the quantum 
spectrum of the underlying system. 

Periodic orbits play the main role in Gutzwiller’s non-integrable systems 
quantization method. The hnal aim of the method consists in the evaluation 
of the density of levels 

n 

in the terms of solutions of classical equations of motion. 

Using the expression 


6{E - E^) 


-^>hm--- 

TT e-s-0 E„ — E — le 


■} 


we could obtain 


PiE) 



E-H 


Operator under the trace is the Green function: 


(4.5) 


G(qA,qB,E) 


E 


E — En 




1 

E-H 




= / S{qA - q)S{qB - q)dq 


And therefore 


G{qA,qB,E) = {qA\- - 

bj — H 


Thus (4.5) could be rewritten in the form 

p{E) = --^SpG 

71 


The further procedure implies the construction of Green’s function semi- 
classical approximation (and then its Fourier transform) and calculation of 
the trace. Gutzwiller had shown that this procedure results in the following 
expression for the levels density: 


OO 

p{E) ^ p{E) + EE Ap^k{E) cos 

p k=l 


k 




(4.6) 
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Here p{E) is the smoothed density of levels that could be obtained via 
Thomas-Fermi approximation or Weyl’s formula for billiards. The sum 
marked by index p is evaluated over all ’’primitive” periodic orbits and sum 
in k — over fc-reiteration of these orbits. The phase for every periodic orbit 
consists of the action along this orbit Ip{E) and Maslov index pp. Amplitude 
Ap^k is determined by the expression 


^p,k{E) 


Tp(E) 

Tik\J det(M — /) 


where Tp{E) = 81/dE is the orbit’s period and M is the monodromy ma¬ 
trix, that is well known from the classical analysis of motion stability. The 
formula (4.6) is called the Gutzviller trace formulc{^ and expresses density 
of quantum spectrum through values that are calculable in the context of 
classical mechanics. At the same time this expression could be understood 
as universal semi-classical quantization condition that is correct both for in- 
tegrable and non-integrable systems: highly excited (semi-classical) energy 
levels are the points where right hand side of the trace formula has poles. 

Although many important results were obtained with the trace formula, 
not all its analytical features are clear enough for now. Mainly this is due 
to the difficulties in the corresponding classical calculations. First of all we 
encounter the problem of periodic orbit evaluation: their number grows ex¬ 
ponentially with growth of period and all of them are by dehnition unstable. 
At the same time there is a problem of adequate description of their contribu¬ 
tions — the summation problem. And hnally the generalization of Gutzwiller 
formula on the considered multi-well case is nontrivial. Thus the problem 
of numerical integration of the Schrodinger equation becomes the basic one 
for calculation of the semi-classical part of the spectrum of quantum systems 
that are chaotic in classical limit. In the following we will briefly review the 
main types of calculation problem arising in research of the manifestations 
of quantum chaos in specihc physical models: evaluation of the energy spec¬ 
trum and investigation of its features, obtaining stationary wave functions 
and their analysis in different representations, modeling of time evolution of 
time-dependent states. We will analyze in detail the main numerical methods 
for these problems — the matrix diagonalization technique and the spectral 
method. 


^Feynman called it one of the main achievements of theoretical physics of the twentieth 
century 
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Chapter 5 

Numerical Methods In 
Multi-Well Potentials 


Studies of deterministic chaos, both classical and quantum, more than other 
domains of modern physics derive their development from computational 
power growth. The number of scientihc papers on that topic published per 
year grows as N ^ where A ~ 0.23 year“^. That value of A is signihcantly 
greater than the growth rate for the total volume of scientific publications 
(A ~ 0.046 year“^), but is very close to the growth rate of global computa¬ 
tional power A ~ 0.25 year“^. And it is not surprising at all because the main 
body of those papers is devoted to quantum chaos researches in numerical 
experiments. 

In the present chapter we develop numerical methods for analysis of quan¬ 
tum chaos problems and demonstrate the advantages of the spectral method 
(SM) in comparison with the matrix diagonalization technique (MD) in ap¬ 
plication to the solution of Schrodinger equation in smooth potential systems, 
in particular with multiple well PPS. 


5.1 The Matrix Diagonalization Technique 

5.1.1 General theory 

Let us consider the Schrodinger equation for a system with discrete energy 
spectrum 

= En^n, UEN (5.1) 
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and let there be full ortlionormal basis of functions 


ipk,k &K, {(pk>\y:>k") — ^k'k" 
where N and IK are countable sets. 

The basis functions ^pk are solutions of another Schrodinger equation 

kp^n ^nPri'/ ^ ^ IK, 


they are given analytically or are obtained numerically in an independent 
way. 

Obviously there exists a decomposition 


a^k^ = (Pkl^n)- 


keK 


The solution of the Schrodinger equation (5.1) by the matrix diagonal- 
ization technique implies the following: 


1. the set IK is presented as a direct sum of the subsets 


K = K © K', 


such as IK is hnite and IK' is a countable set. 


2. original Hamiltonian of the problem (5.1) is presented in the form 


H = H + H', 


where by dehnition 

{ipk'\H\(pk'') = {(pk'\H\ipk"), k',k” e IK (5.2) 

and all other matrix elements of H are zeros. 

3. the eigenvalue problem 

= En^n, ^ € K (5.3) 

is solved numerically, where 

a^k^ = (Pkl^n), (’hn'I'^n") = W'- 

fceiK 


(5.4) 
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Let us find the conditions under which the numerically obtained En and 
On are good approximations to the original En = En + E'^ and d'n = + 

or in other words the conditions for smallness of E'^ and 
Let us rewrite ( |5.1[ ) in the form 

(H + H') (f „ + <I>y = (B„ + E'^) (■!>„ + 'ry 

and simplify it using ( |5.3[ ) to obtain 

H'^n - E'J>n + HK - En% + H'^'^ - E'Jl'^ = 0. (5.5) 


Taking scalar product of d/^, in (5.5), and taking into account that, according 
to (5.2) and (5.4), 

{^n\H'\^n)=0, |'hnr = l, 

we obtain 


{^n\ 

nO CD TT/ 

/ LLk/ tLuf/llyyf _ 

_ k'eK _ (TJ 

H'\ 

^n) 

\ 

1 ^ n 1 


E ay'al?’ 1 - 

IT' 

1 n 

IV 2 


fc'eK 


It is easy now to see, that smallness of 

4"^ = keK' (5.6) 


is a sufficient but not necessary condition for smallness of E'^. In order to 
hnd the conditions for smallness of (5.6) it is convenient to make use of the 
Wigner representation 


q) = I dxe-^nP-^n {q + f) {q - f) 
q) = I dxe--^P^ipk [q + I) (Pfc (g - f) , 


and then we have 

^ o [ dpdq¥^\p,q)(p^^\p,q), 

{27in) J 

where D is the dimension of the conhguration space. 

According to the principle of uniform semiclassical condensation of quan¬ 
tum states |12] in semiclassical limit, the Wigner functions d/i^^(p, g) and 
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P^\piQ) ci-re localized on energy surfaces H{p,q) = and h{p,q) = 
respectively 

{P, q) ^ ^ {H{p, q) - En) , {p, g) ~ 5 {h{p, q) - Ck). 


The coefficients ( 5.6[ ) are negligibly small when the corresponding energy 
surfaces do not have common points, and therefore analysis of applicability 
of the matrix diagonalization method is reduced to analysis of the classical 
phase space of the system under consideration |13]. 

5.1.2 ID harmonic oscillator 

As the simplest example let us consider calculation of the energy spectrum 
of a one-dimensional harmonic oscillator 


H{p,q) = ^ + ^^ 

on the basis of eigenfunctions of another harmonic oscillator 


(5.8) 


hr \ P I 2r 

h{p,q;u) = — +u 

The basis functions, corresponding to energy levels 


have the well-known form 


Ck = hu I k + - 


where 


/u \ i ffkiO -4 . 


(5.9) 


(5.10) 




[n/2] 


^n(0 = (-l)V^;^e-«^ = 5: 


~ 2m)! 


(20 


n—2m 


are Hermit polynomials and [n/2] means integer part. 


The Wigner form for the basis functions (5.10) reads 

H-OO 


Pk{p,q;uj)= / dxe (g + |;w) w) = (-l)^2e 
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where 


^2 ^ 


fun 


and Lk are Laguerre polynomials. 


Taking a; = 1 in the expressions for energy levels (5.9) and basis fnnctions 


(5.10) we immediately obtain the exact resnlt for the original problem 


En = h ( n + - 


bnt for arbitrary u diagonalization of the original Hamiltonian (5.8) in the 


chosen basis represents already non-trivial procednre, becanse the corre¬ 
sponding matrix is non-diagonal 


Hmn 2 a;) 2 ) 

{n -|- l)(n -|- 2)Sm,n+2 + \/{m + l){m -|- 2)6n,m+2 


As the original and the basis fnnctions are semiclassically localized in the 
phase space on the ellipses \ \ = En and \ respectively, 

simple geometrical analysis of the intersection conditions for the correspond¬ 
ing manifolds dehnes the basis set, which is necessary and, as the nnmerical 


experiment conhrms (£g.5.1), snfficient, to correctly determine the reqnired 
state En 

uj, — \ En < h \ k -\— ) < max (ca, — ) En- 


mm 


bJ 


Ul 


As can be seen on Fig,5.1, at hxed set of basis fnnctions ki < k < k 2 


the region of compntationally stable resnlts for energy levels represents the 
interior of a cnrvilinear tetragon formed by a pair of hyperbolas and a pair 
of straight lines: 


hu (ki + I 


< En < hw [k2 + I) 
<En<^{k2 + i) 


(5.11) 


It is convenient to introdnce new variables 


7] = \nE, ^ = Inw, z/ 1,2 = In ( fci ,2 + - 


(Figg 


in which the stability region (5.11) transforms into the interior of the sqnare 


Ul < T] ± ^ < 1/2. 


( 5 . 12 ) 
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a b 




Figure 5.1: Dependence of numerically found energy levels of harmonic os¬ 
cillator (5.8) on the auxiliary basis frequency oj for 0.03 < uj < 10.0, ki = 2, 
^2 = 100 (a) and 0.75 < u < 1.25, ki = 400, k 2 = 500 (b): points represent 
the numerical data, thin solid lines are exact energy levels and thick solid 
lines limit the computational stability region determined from semiclassical 
analysis (5.11) 
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a b 




Figure 5.2: Dependence of numerically found energy levels of harmonic os¬ 
cillator (5.8) on the auxiliary basis frequency oj for ki = 2, ^2 = 60 (a) and 
ki = 130, k 2 = 190 (b): points represent the numerical data, solid lines are 
exact energy levels and gray color fills the stability region determined from 
semiclassical analysis (5.12) 
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Obviously the optimal basis will be that which corresponds to maximum 
vertical dimensions of the stability region. In the considered example it 
evidently corresponds to frequency a; = 1. 

An even simpler basis set can be built on plane waves — eigenfunctions 
of one-dimensional billiard, i.e. inhnitely deep potential well with walls in 
the points Xi ^2 = ±a: 


(Pk{x;a) = 


±sm^{x + a), \x\<a 


0 , 


\x\ > a 


(5.13) 


The semiclassical phase space localization analysis can be applied to the basis 


functions (5.13) as well, if we take into account that they are localized on 


the straight lines 


X = ±a 

P=^fak 


We should note that, unlike the harmonic oscillator, the plane waves basis 


(5.13) cannot be cut from below: ki = 1. As result the stability region on the 


(£^, a)-diagram for computation of the harmonic oscillator spectrum (5.9) in 


the plain waves basis (5.13) with indexes A; = 1,..., ^2 is represented by the 


interior of a curvilinear triangle with base E = 0 and bordered by a parabola 
and a quadratic hyperbola, depending on k 2 . 


En<i 

En < 


8a2 


In the variables 

(h = log 2 E, i = log 2 a, u = log 2 Trhk 2 ) 
the computational stability region takes the form of a normal triangle (see 


£g,5.3r 


7 ] <2^-1 

I) < 2(1/ - 0 - 3 


(5.14) 


The vertex of the triangle (5.14[) 


e = 


u-l 


,r] = u 
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corresponds to the optimal basis with parameter 


a = 


Tihko 


(6.15) 


which allows correct calculation of energy levels up to 

F — 1? 

■^max . 


Therefore the plane waves basis ( |5.13 ) with indexes A; = 1,..., ^2 allows the 
correct calculation of the harmonic oscillator energy levels (5.9) with indices 
= 0,..., A ;2 = iV (see Figj5.3|3), where 


n 


N ='^k2--. 

4 2 

For large basis dimensions the fraction of correctly calculated energy levels 
tends to the limit 

, 77 TT ^ , 

lim — = - = 75% 5.16 

fc2^oo k 4 


5.1.3 Quartic oscillator in ID and in 2D 


As a more complicated example let us consider calculation of energy levels 
for the quartic oscillator 

H{p,x) = + p‘^. (5-17) 

As distinct from the harmonic oscillator, this quantum mechanical prob¬ 
lem does not allow analytical solution, but there is a very accurate semiclas- 
sical approximation for the quartic oscillator spectrum [H2] 


TT 


2 f 3^/2(n+|) 


5 I 11 
8"^ 6 


'(i)‘ 


97r(n+|) 




nij 


+ 


1 + 


ill 

Ttt 


97r(n+|) 


(5.18) 


97r(n+i) 




Formula (5.18) gives the energy levels with 15 correct digits, which is quite 
sufficient to check the numerical results in the given case. The results of 
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parameter a: points represent the numerical data, solid lines are exact energy 
levels and gray color fills the stability region determined from the semiclassi- 
cal analysis (5.14) b) Dependence of absolute error on main quantum number 
n for optimal basis parameter (5.15). The vertical line corresponds to the 
applicability limit (5.16) of the matrix diagonalization method. 
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Figure 5.4: Calculation of energy levels of quartic oscillator (5.17) in the 
plane waves basis (5.13) for ki = 2, k 2 = 2500. a) Dependence of numeri¬ 
cally found energy levels on the auxiliary basis parameter a: points represent 
the numerical data, solid lines are exact energy levels and gray color hlls the 
stability region determined from the semiclassical analysis b) Dependence 
of absolute errors on main quantum number n for optimal basis parame¬ 
ter. The vertical line corresponds to the applicability limit of the matrix 
diagonalization method. 
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Figure 5.5: Calculation of energy levels of quartic oscillator (5.17) in the 
harmonic oscillator basis (5.10) for ki = 2, ^2 = 250. (a) Dependence of 
numerically found energy levels on the auxiliary basis frequency uj: points 
represent the numerical data, solid lines are exact energy levels and gray 
color fills the stability region determined from the semiclassical analysis b) 
Dependence of absolute errors on main quantum number n for optimal basis 
parameter. 
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such a check are presented on Fig j5.4| for the plane waves basis (5.13) and on 
Fig j5.5| for the harmonic oscillator basis ( |5.10 ). 

Semiclassical phase space analysis gives the following results for the quar- 


tic oscillator (5.17): the stability region is determined by the conditions 


En < 

E < 

^ 4^2 <^2 


for the plane waves basis (5.13) and 


En<2nw{k + \)-'^ 

En<S(k+C 


in the harmonic oscillator basis (5.10). Now it is easy to obtain expressions 


for optimal basis parameters and maximum energies of correctly calculable 
states: 

„ ( nhk \ 3 rpPW ( nhk \ 3 

® ’ ^max - 

w^[2h{k + yi, E«y = i[2h{k + yK 

Taking into account ( |5.18[ ), we obtain for the relative number of correctly 
calculated energy levels 


N 

k2 


£!(i) 

6^^ 


0.874 


for the plane waves basis (5.13) and 


N 

k2 


Elli) 

27r533 


0.897 


in the harmonic oscillator basis (5.10). Hence it follows that the latter is 
slightly more efficient. 


Harmonic oscillator (5.10) and inhnitely deep potential well (5.13) in fact 


exhaust the set of exactly solvable one-dimensional quantum systems whose 
eigenfunctions can be used as a basis for matrix diagonalization. There are 
many more possibilities in the models with dimensionality of more then one. 
Further, for simplicity we consider only two-dimensional systems, but the 
results can be trivially generalized for higher dimensions. 
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The simplest type of two-dimensional basis can be constructed from prod¬ 
ucts of eigenfunctions of exactly solvable one-dimensional problems, for ex¬ 
ample, from plane waves 


Pkx,ky[X,y, Clx, Uy) 


^JYxOty 


: sm 


2ar 


[x ax) sin 


nky 

2a,, 


[x + a. 


y)i 


(5.19) 


harmonic oscillator eigenfunctions 


Pkx,ky (^7 ^X7 ^y) 


yOJx^y 


,iHh 




^-(ujxx'^+ujyy'^)f{2h) 


y/nh kx\ky\ 

or as a combination of both of them 


(5.20) 


Pkx,ky 


{x,y,a,u) = ^ 


sm 


irkx 

2a 


{x + a) 


yJVMe-^ 

y/2.Gky\ 


(5.21) 


Efficiency of the basis types ( 5.19[5.20[5.21 ) depends on the problem under 
consideration, but in any case simple semiclassical phase space analysis al¬ 
lows us to choose the most preferable among them. Skipping cumbersome 
exact expressions, let us give only the ultimate results (Fig j5.6[5.7^ ) of the 
semiclassical optimization of the spectrum calculation in the two-dimensional 
potential of coupled quartic oscillators (CQO) 


U 4 {x, y, a) = x'^ + y"^ + ax‘^y‘^. 


(5.22) 


Quite often the Hamiltonian of the system under consideration has a dis¬ 
crete symmetry. For example, the CQO potential (5.22) is invariant under 
transformations of square symmetry group In such cases it is conve¬ 

nient for many reasons to calculate the states of different symmetry types 
independently. Firstly, for investigation of statistical properties of energy 
spectra in quantum chaology we have anyway to exclude pure spectral series 
— sequences of states with one and the same symmetry type. Secondly, as 
a rule, states of different symmetry, even and odd for example, are usually 
very close in energy if not degenerate. Therefore numerical computation of 
all symmetry types of states together leads to very ill-conditioned matri¬ 
ces, while exclusion of certain symmetry type improves the conditionality. 
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Figure 5.6: Calculation of energy levels of CQO potential (5.22) with a = 6 
in the harmonic oscillator basis (5.20). a) Dependence of numerically found 
energy levels on the auxiliary basis frequency u = = ujy'- points represent 

the numerical data, solid lines are exact energy levels and gray color fills the 
stability region determined from the semiclassical analysis b) Dependence of 
absolute errors on main quantum number n for optimal basis parameter. 
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And lastly, computation of different symmetry types one-by-one runs evi¬ 
dently faster than determination of all the states altogether. For example, 
the basis for determination of Ai-type states in the CQO potential (5.22) is 
constructed from symmetrized combinations of the form iprnn + y^nm of basis 
vectors (5.19) with Oj, = Uy, {m,n odd) or (5.20) with = ojy {m,n even). 

The characteristic feature of polynomial potential is the sparse band 
structure of the Hamiltonian matrix (Fig j5.7[ 3) in the basis of harmonic os¬ 
cillator (5.20). For example, for large basis dimensions n the bandwidth for 
Ai-type states in the CQO potential ( |5.22| ) is equal m 2-\/^. Clearly for 
polynomial potentials it makes sense to use special routines for band ma¬ 
trix diagonalization which allows us to economize both the CPU time and 
memory usage for computations. 

Simple analysis shows that the number of non-zero matrix elements for 
Ai-type states in the CQO potential (5.22) never exceeds lln. It means that 
the used basis vectors ordering is not optimal and that there possibly exists 
another ordering which leads to band matrices with constant bandwidth m = 
5, or at least with slower growing bandwidth. Such ideal ordering would 
considerably speed up the computations but its search represents a non¬ 
trivial task. 


5.2 The Spectral Method 

The spectral method (SM) for the solution of the Schrodinger equation was 
proposed in the paper m in application to ID and 2D potential systems, 
but it can be easily generalized for the Schrodinger equations of arbitrary 
dimensions: 


+ U{xi,...,xd) 

^=1 


'Ipni.Xi, . . . , Xd) = En'llJni.Xi, ...,Xd), 


where D is the dimensionality of the system conhguration space. Let us as¬ 
sume that the potential U{xi,... ,X£)) allows only hnite motion for all ener¬ 
gies, therefore our task is to hnd discrete energy spectrum En and stationary 
wave functions 'ipn{xi,... ,xd)- 

Let us consider time-dependent solution 'ifj^xi,... ,XD]t) for the corre- 
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Figure 5.7: a) Comparison of numerically found energy levels of CQO poten¬ 
tial (5.22) with a = 6 in the plane waves basis (5.19) and in the harmonic 
oscillator basis (5.20) for optimal parameter of both bases, b) Distribution of 
non-zero Hamiltonian matrix elements in the harmonic oscillator basis (5.20) 
for Hi-type states in the CQO potential (5.22). 
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spending non-st at ionary Schrodinger equation 

with some in principle arbitrary initial condition 


tlj{xi, ...,XD]t)= ihdtiljiyXi, ...,xdG) 


ijQ^Xi, Xd) = 'ipni.Xi, ...,XD]t = D) 

. Applying the decomposition 

OO 

i)o{Xi, . . .,Xd) = y^^ani>n{xi, . . .,Xd), 

n=l 


we obtain 

OO 

^jJ{Xl, ...,XdG) = y^^an'tpnjxi, ■ ■ .,XD)e~'~^. 

n=l 

Here and further we imply that the wave functions 'ipnixi,... ,xd) are or¬ 
thonormal 


j dxi... dxoi’iixi ,..., XD)iJk{xi ,..., Xd) = dik- 


Let us consider autocorrelation function of the form 

/ OO 

dxi... dxD^oixu • • •, XD)i’{xi ,..., Xd; t) = ^ |a„|^ e"*" 

n=l 

We assume the initial wave function 'ipQ^xi, ... ,xd) to be normalized 

p OO 

f dxi . . .dxD'lpoiXi, . . . ,Xd)'4’o{Xi, . . . ,Xd) = ^ \anf = 1, 


(5.23) 


n=l 


so that P{0) = 1. The Fourier transform of the autocorrelator (5.23) contains 
information about the energy spectrum En of the system 
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For determination of stationary wave functions 'ipn{xi,... ,xd) we will need 
the Fourier transform of ... ,XD',t) itself 


^jJ{Xl,...,XD;E)= / dtc ^p{Xl,. .. ,XD;t) = h'S^an'ijJnixi,... ,Xd)HE - En) 


n=l 


(5.25) 

Naturally in practice we never try to find all En and tl’nixi,... ,xd)- Usu¬ 
ally the task is to determine all the energy levels En inside a given interval 
El < E < E 2 with certain accuracy 5e and to calculate the correspond¬ 
ing stationary wave functions 'ipn{xi,... ,xd) on some finite set of points 


-|- kiAxi, ki = 0,..., Ni> also with hnite accuracy. Further we 


Ck,) _ .^(0) 

will assume for simplicity that all Ni are equal. 

Let us show how to apply expressions ( 5.24[5.25 ) for computation proce¬ 
dure construction. It is evident that in practical calculations we always deal 
with only a hnite number of data known also with hnite accuracy. In our case 
it means that the principal function — ‘ipixi,... ,X£)]t) — will be known on 
a hnite set of points both on temporal and spacial coordinates. Accordingly 


the autocorrelation function (5.23) can be calculated also only on a hnite set 


of points tk = kAt, k = 1,..., M, and therefore we have to apply the hnite 


analogues of the Fourier transforms (5.24) and (5.25) 


Pt(E) = y 


/ P{t) = 

i n=l 


'AriE-En) 


(5.26) 


i>T{xi,... ,xd-,E) = ^ J dtc’^^ 'ijjixi, ...,XD]t) 

0 

CXD 

= an'lpnixi, ..., XD)dT {E - En) , 


n=l 


(6.27) 


where T = MAt and the hnite analogue of h-function takes the form 


St{E) = y f dtc’n 

Mx) = 

As distinct from the usual 5-function, 5r(0) = 1 (see Fig j5.8^ ). Therefore 
each energy level En corresponds to a sufficiently sharp peak — local maxi¬ 
mum of the function \Pt{E)\, situated at the point E = En (Fig. |5.8|3). The 


iET 

e TT-l 
n 


fr tD 
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Figure 5.8: a) \fT{x)\ profile (dashed line) and \fH{x)\ pro¬ 

file (solid line). b) \Pt{E)\ (dashed line) and \Ph{E)\ (solid 
line) for the model system with energy spectrum F'i, 2,...,9 = 

{—30, —10, —7,10,11, 20,20.7,40,40.63} (marked by solid vertical lines) and 
ai, 2,...,9 = {0.5,0.5,0.5,0.25,0.25,0.2,0.2,0.15,0.15} with h = T = 1. 
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typical profile of the function \Pt{E)\ has many local maxima (Fig j5.9[ ), but 
a major part of them has absolutely nothing to do with the energy levels of 
the system under consideration — they appear due to the oscillating charac¬ 
ter of 5t{E). So the full number of maxima on Fig j5.8[ 3 equals 94, while the 
corresponding system has only 9 physical energy levels in reality. Therefore, 
formally analyzing \Pt{E)\, we will obtain plenty of extra ’’parasite levels”. 
The unique characteristic that allows us to distinguish such phantom levels 
from the real ones is the relative smallness of \Pt{E)\ in corresponding local 
maxima. However if we were just to ignore all the levels that have \Pt{E)\ 
peak amplitude less than a certain hxed threshold value, we risk loosing 
some real physical levels, which correspond to small values of |an|. In this 
connection it is very useful to apply the weighted Fourier transform: 


P^{E) = ^f dtP^P{t)w{t) = E l«n|'E [E - En) 

0 n=l 

T 

‘ipwixi,... ,Xd]E) = ^ J dtpxr'ijjiyci,... ,Xd] t)w{t) = 

0 

oo 

= E an^/jnixi, . . . , {E - En) 

n=l 

Su,{E) = y f dtP^w{t) 

0 


(5.28) 


where the weight function w{t) satishes the conditions 


ta(0) = w{T) 

T 

J w{t)dt = 1. 


(5.29) 


Therefore E(0) = 1 for any such w{t). The simplest function, satisfying the 
conditions (5.29), is Hann function 

2TTt 

WhR) = 1 - cos—. 


for which 


dH{E) = 


iEt 

eXT-l 


“T ,Et o_ 


• Et I o 


fni 


_ sin(7r3;)e^^^ _ frjx) 

— 1—X‘^ 


= fx (f) 
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E 


Figure 5.9: \Pt{E) \ for calculation of i7-type energy levels in the quadrupole 
oscillation potential. 
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Inclusion of the weight function in the modihed Fourier transform al¬ 
lows us to diminish the relative amplitude of the phantom peaks in \Ph{E)\ 
(Fig jS.S^ ) and even slightly decrease their number — to 76 for \Ph{E)\ in¬ 
stead of 86 for \Pt{E)\ (Fig j5.8| 3). Numerical analysis of the |/r(a;)| shape 
shows that the second maximum is situated in points XTmax ~ ±1.43 and 
has amplitude frmax ~ 0.217, while for the |//f(a;)| the analogous anal¬ 
ysis maximum is situated almost twice farther from the principal one — 
XHmax ~ ±2.36, and its amplitude appears less for whole order of magnitude 
— frmax ~ 0.0267. Amplitude of more distant maxima decrease even faster 
(FigjK^). 

Besides that, an important independent criterion to estimate the accuracy 
of numerical results can be realized using the semiclassical approximation for 
the stair-case states number function n{E): comparison of numerical n{E) 
with the semiclassical one allows us to roughly estimate the number of lost 


or acquired extra levels (£g,5.10). 

Analyzing the positions and amplitudes of the \Ph{E)\ local maxima one 
can determine the energy levels En and absolute magnitudes of the coeffi¬ 
cients On 

l^nl = \/\^H{En)\- 


If in the expression 


...,xd; En) = ttnfJnixi ,..., xd ) ± ^ • • • , XD) 5 H{En - Ek) 

k^n 

we now neglect the terms containing SniEn — EQ (compare with 
can determine the eigenfunctions up to a phase factor: 

'ipnixi, ...,xd)^ ...,xd; En) (5.30) 

|®n| 

Here and forth we consider that there are no degenerate levels among En¬ 
in practice in any system with a degenerate spectrum we can (and should) 
remove the degeneracy by a certain choice of the initial state ffo^xi,... ,xd)- 
However such an approach is applicable only while the given energy level 
En is sufficiently separated from the neighbors. Indeed, if one or more levels 
are situated at too short a distance, the position of corresponding maxima 
\PH{En)\ considerably differs from the actual values En- Moreover for suffi¬ 
ciently close levels the {PniEQl prohle will have only a single common peak 
(Figj^). 


5.28), we 
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E 


Figure 5.10: State number n{E) (dashed line) in Thomas-Fermi approxi¬ 
mation versus the numerically obtained n{E) (solid line) for lower umbilic 
catastrophe potential. 
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El 

ai 

Exact 


-30.00 

0.50 

From 

\Pt(E)\ 

-30.03 

0.50 

From 

\Ph{E)\ 

-30.00 

0.50 



E^ 

04 

Exact 


10.00 

0.25 

From 

\Pt{E)\ 

9.89 

0.25 

From 

\Ph(E)\ 

9.53 

0.21 



Ej 

Idyl 

Exact 


20.70 

0.20 

From 

\Pt(E)\ 

21.10 

0.19 

From 

\Pm(E)\ 

- 

- 


E 2 

I 02 I 

E 3 

hal 

-10.00 

0.50 

-7.00 

0.50 

-10.09 

0.50 

-6.91 

0.50 

-9.97 

0.50 

-7.03 

0.50 

E, 


E, 

lael 

11.00 

0.25 

20.00 

0.20 

11.24 

0.27 

20.22 

0.19 

11.47 

0.21 

20.35 

0.18 

Es 


Eg 

lagi 

40.00 

0.15 

40.63 

0.15 

40.32 

0.16 

- 

- 

40.31 

0.15 

_ 

_ 


Table 5.1: Energy spectrum 
system from Fig. 5.8)3 


and eigenstates amplitudes a„ for the model 


In general, it looks impossible to determine exactly the minimum sepa¬ 
ration between close levels at which they still can be resolved, but we can 
state with conhdence that this value is of the order of AEmin = h/T — the 
natural width of level. Therefore all duplets and multiplets with separation 
less than AEmin will give only single local maxima, and some real levels will 
be lost in the computed spectrum (see Table 5.1). 

From Table 5.1 we can see that analysis of \Ph{E)\ gives considerably 
more accurate results for isolated levels. On the other hand, \Pt{E)\ has a 
smaller natural level width and gives more adequate results for close levels. 
However the preference of \Pt{E)\ in resolution of close levels is so insignih- 
cant that \Ph{E)\ appears more preferable in the majority of cases. 

Taking into account that P{t) as well as ... ,XD]t) is calculated 

only on hnite set of points tk = kAt, k = 0,..., M, the most natural way to 
obtain Ph{E) is the discrete Fourier transform 


^ IVl A , IVl 

PniEk) = -^^P{tj)w{tj)e "m = — '^p(^tj)w{tj)P^, (5.31) 

j=0 j=0 

where Ek = yk, k = —M/2,..., M/2. Therefore Ph{E) appears to be 
calculated in energy interval — ^ < E < ^ with step AE = h/T exactly 
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equal to natural level width. 

As w{T) = t(;(0) = 0, the discrete Fourier transform (5.31) in fact coin¬ 
cides with the formula for numerical integration by the trapezium rule|45j. 
applied to (5.28), and it is easy to estimate its error 


RniE) = 


TAB 

12 




t=T 


d2 

dP 


where 0 < r < T. But 

^ |a^|2 7g27ri|-(£!-E„) _ _ 1^2ni^(E-E„-{^)\ 

- ^ [E — En + y) — 7^ yj-j — J-Jn — y) ^ 

< P E K? {E - Er,f + \{E-E^ + hY + \{E-E^-hY 

+ (f )1 = ^ + 2 (^) V (f )^ 


< 


hyp ^—2m^ 


n=l 


< 


n=l 

oo 

~ E 

n=l 


where E and De are respectively mean energy and dispersion in the initial 
state 'ipoi^i,... ,xd)- 



E = J dxi... dxD'ipoixi, xd)H' 4 !o{xi, xd) 
Jdxi... dxoi’oixi, ...,xd){H- Ef' 4 )o{xi, ...,xd) 

H = -^f:df + U{x,,...,XD), 


= E Kl^Er, 


n=l 

oo 


n=l 


and we hnally get 


TAR p - 2 2 1 e'^AR 

Rh{E) < ^j _2 \{E — E) -|- Y 


Qh? 


3T 


Therefore formally integration in (5.28) allows us to calculate Ph{E) for 
any energy values, but in fact applicability of such an approach is limited by 
the energy region 
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which is definitely not better than for (5.31). Application of the approxi¬ 
mated integration formulae of higher orders also does not give any improve¬ 
ment to (5.31), because the error estimate for numerical integration of (5.28) 
by n-th order method reads 


<>(£) 


TAD 

nlh'^ 


(E - Ey 


For calculation of , xd'R) at tk = kAt, k = 1,...,M we can 

apply step-by-step the split operator method 


ip{xi,... ,XD]t+At) = e 


h At q2 


e »=i 


y{xi,... ,XD;t)+0{AC). 


0 D 

iAAI y, df 

Action of the differential operator e ” is also calculated with the 

help of the discrete Fourier transform 


y D 

e ...,XD]t) = 


Ni/2 

= E 


N 2/2 Nn/2 

E ••• E 


y. 


nin 2 ...nr> 


(f)e*= 


° r 

E 




1 ) 


ni=—Wi/2 772=—A^2/2 n£)=—N£}/2 

Ni/2 N 2/2 Nn/2 _2m Y 

ynin2...nn{t) = -D^ E E ••• E y{xi, . . . ,XD]t)e 

n ni=—Vi /2 n2=—V2/2 nn=—Nn /2 

fc=i 


With sufficiently large numbers of time steps M and coordinate grid nodes 
N the results of computations by the spectral method really do not depend 
on the arbitrary initial wave function yo^Xi,... ,xd), but it is a reasonable 
choice of such a function that is the most powerful means for the optimization 
of computations by the spectral method in order to achieve sufficient accuracy 
of the results with minimal expenses of computational power. 

It is convenient to choose the initial wave function ^oixi,... ,X£,) as a 
linear combination of Gaussian wave packets of minimal uncertainty 


yoixi ,.. .,xd) = e S 




+i 


Piixj-XiQ) 


(5.32) 


Combining functions of the form (5.32) with different parameters a, Pi and Xio 

,xo) with the desired 
Its position 


it is possible to selectively excite eigenstates i/Juixi, 
symmetry properties and lying in the desired energy interval 


{E) = {yoixi, ■ ■ . ,Xd)\H\'4!o{xi, . . .,Xd)) 
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and approximate width 


.. .,xd)\{{H - {E)f\^jJo{xl,.. .,xd)) 


are determined by parameters of the Hamiltonian and Gaussian wave packets 
that generate the initial state 


'Ipoixi, . . .,Xd) = y^^an'lpnjxi, . . .,Xd), 

n 


and can be varied in a wide range if desired. This in particular allows us to 
calculate low-energy and high-energy states separately, which leads to more 
efficient distribution of the computational efforts. 


In the spectral method computations the preponderant fraction of CPU 
time is invested in the calculation of '4^{xi,... ,XD]t), and namely to multiple 
direct and inverse discrete Fourier transforms. For two-dimensional prob¬ 
lems the characteristic calculation time-scales are MN'^lnN which is sub¬ 
stantially better than for the matrix diagonalization method. Another 
important advantage of the spectral method is the possibility to calculate 
multiple eigenfunctions in parallel computations by (5.30), which allows us 
signihcantly to economize the CPU time. For example, the computation time 
for 100 eigenfunctions is only twice longer that for one eigenfunction. 


The unlimited possibilities of the spectral method in computational ac¬ 
curacy improvement are demonstrated in Fig, 5.11 for a problem of determi¬ 
nation of two close energy levels in the quadrupole potential. With time step 
number M increasing from 2® to 2^®, peaks of P{E) became more and more 
pronounced (Fig jS.ll^ ). For M = 2® and M = 2^° energy resolution is yet 
too small, and two neighboring levels look as one. At M = 2^^ the duplet is 
already resolved but obtained values for energy levels differ signihcantly from 
the real ones. Sufficient accuracy of energy level determination is achieved 
at M = 2^^ and further increasing of M does not lead to remarkable changes 
in calculated energy levels values (Fig jS.llj p). The accuracy of obtained re¬ 
sults grows very fast with increasing of M and can be made arbitrarily high 


(Fig j5.11| 3). 










5.3. COMPARATIVE ANALYSIS OF MATRIX DIAGONALIZATION AND SPECTRAL METHODS 


(a) E (b) log^AE (c) 
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Figure 5.11; Determination of two close energy levels En and En+i in the 
quadrupole oscillations potential with hF = 18 and different time step num¬ 
bers M = 2^, k = 9,10,..., 18: a) \P^^\E)\ in logarithmic scale; b) values 
of En and En+i, c) consecutive corrections En — En^ and 


5.3 Comparative analysis of matrix diagonal- 
ization and spectral methods 

The method of Hamiltonian diagonalization is the most traditional way for 
numerical solution of the Schrodinger equation. The spectral method for 
solution of the same problem represents in its turn a newer one and for many 
reasons a more preferable approach. 

One of the most fundamental disadvantages of the matrix diagonaliza- 
tion technique is the rather poor choice of exactly solvable models whose 
eigenfunctions can be taken as a basis for subsequent diagonalization of the 
Hamiltonian under consideration. As a rule, properties of the Hamiltonian 
pose very rigid limitations on the auxiliary basis parameters, therefore in 
most cases of matrix diagonalization implementation, only one free parame¬ 
ter remains — it is the auxiliary basis dimension N. Further it is necessary 
only to determine the minimal dimensionality N sufficient for the achieve¬ 
ment of desired resulting accuracy. Such simplicity of the matrix diagonaliza- 
tion method results in its insufficient flexibility: in practice the application 
of matrix diagonalization is justihed only for those potentials that can be 
approximated at least locally by some exactly solvable model. But such a 
limitation cannot be satisfied for many important problems, especially in the 
potentials with many local minima. 

The spectral method uses a natural basis of free particle wave functions 
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— such a basis is equally good, or better to say equally bad, for potentials of 
any form. Such fundamental indifference of the spectral method to shapes of 
potential energy surface is the main reason of its universality. Compared to 
matrix diagonalization, the spectral method has much more flexibility — the 
researcher is free to choose both the length and step of the computation grid 
in time (T and At) as well as in coordinate space {Li and Axi). In the same 
time choice of the nodes number N is limited by the computational efficiency 
requirement: the applicability condition for the fast Fourier transform algo¬ 
rithm — the main basis of the spectral method efficiency — assumes that all 
Ni do not contain large simple factors; ideally all of them should be integer 
powers of two {Ni = 2^*). And, last but not least, the main freedom lies in 
the choice of the initial state for the spectral method computations. As it is 
very difficult to give any general recommendation on that point, the spectral 
method computations have become a real art rather than plain technique, re¬ 
quiring great experience and constant practice. Because the spectral method 
is not standardized up to the present time, the fast Fourier transform rep¬ 
resents the only one ready-to-use ingredient for its realization, available in 
many well-known software libraries. Other stages of compntations reqnire 
rather extended although principally simple software development. On the 
other hand, the spectral method algorithm itself can be easily generalized 
for problems of any dimensions, which is not the case with the matrix diag- 
onalization techniqne — reasonable constrnction of hnite mnlti-dimensional 
basis from one-dimensional eigenstates always represents a non-trivial task 
because of the basis vectors ordering problem. 

A quantitative measnre of the nnmerical method efficiency is the growth 
of computational expense — CPU time and RAM usage — with increase 
in the results both quantity and quality. It is useful to compare the matrix 
diagonalization and the spectral method efficiencies in calculation of n energy 
levels of a qnantnm system with hxed relative error e. Taking into account 
the fact that for studies of statistical properties the calculated spectrum is 
inevitably unfolded, it is reasonable to dehne the accuracy as the maximum 
ratio of absolute error of the computed energy levels to mean level spacing 

e=N = p(e)Se. 


For most potential systems the level density p{E) grows qnite fast with 
the energy: 

p{E) ^ E^^-\ 
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where D is the system dimensionality and A is close to unity (it exactly 
equals unity for a harmonic oscillator). Therefore the condition to achieve 
the desired accuracy will be the most critical for levels with maximum energy, 
while the lowest levels will be obtained with higher accuracy than needed 
— this inconvenient feature is due to the very nature of smooth potential 
systems and is equally shared by both methods under consideration. 

In the matrix diagonalization method the absolute computational error 
for sufficiently low energies does not exceed the round-off errors: 

~ SoEkN^, k < r]N, 0 < r] < 1 

where N is the basis dimensionality, rj represents the problem-dependent rate 
of correctly calculable states and (5o is the machine round-off error ((5o ~ 10“^® 
for standard double accuracy numbers). For such levels the relative error 
e ^ appears to be negligibly small and for required basis dimensions 

we get the condition; 

_ n 
Emd = —■ 

V 

In the general case of matrix diagonalization the computation time scales 
as Tmd ~ = n^/rj^ and memory usage is Mmd ~ = n? /ff- However 

in many important particular cases, for example for polynomial potentials, 
the matrix diagonalization implementation involves band matrices, and in 
such cases Tbmd ~ and Mbmd ~ N. 

In the spectral method the computational accuracy is determined by the 
size of time computational grid: 5e ~ l/T, while the time step determines 
the spectral bandwidth where the levels can be determined; At ^ 1/Emax- 
Therefore the maximum relative error in the spectral method scales as e ~ 
E^^-^/T, and the required time step number Nt = T/At ^ n/e. Because 
the spectral method actually uses the plane waves decomposition for the 
computed wave functions, the required nodes number for the computational 
grid for each of the dimensions, equal to the effective basis vector number 
on the same degree of freedom, is determined in the same way as for the 
corresponding basis in the matrix diagonalization method: 



where rjpw is the rate of correctly calculable levels on the plane waves basis. 
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Therefore the computation time scales as 


Tsm-NtN^ In N, 



In 


n 


eripw Vpw 


and memory usage scales as 


Msm ^ N[’ ^ - 

hPW 

As a result the spectral method is generally more efficient for both CPU time 
and RAM usage criteria: 

Msm rj^ Tsm , h 

Mmd nhPw' Tmd p^Vpw Vpw 

The last but not least advantage of the spectral method lies in the fact 
that numerical simulation of wave packets temporal dynamics is included 
into the method algorithm as an auxiliary procedure. In some important 
applications this simulation represents the ultimate goal of research, while 
the determination of energy levels and stationary wave functions is not re¬ 
quired. In such cases an important advantage of the spectral method is the 
possibility of achieving the goal by the shortest way, while it would require 
much more computational efforts to reproduce the same results by the matrix 
diagonalization technique. 



Chapter 6 

Statistical Properties of Energy 
Spectra 


The energy spectra represent historically the first object for investigation of 
quantum signatures of chaos in Hamiltonian systems. Substantial progress 
in the detection of quantum manifestations of classical stochasticity in the 
1980s was connected with the transition to studies of statistical properties of 
energy spectra. Up to that time there had been established the connection 
between spectral properties of complex systems (for example, atomic nuclei) 
and ensembles of random matrices, which created the background for the 
understanding of the complicated pseudo-random nature of energy spectra. 
This in turn stimulated the transition from investigation of the behavior 
of separate levels to statistical characteristics of the energy spectrum as a 
whole. In the research process it became clear that only local (not averaged) 
statistical properties of energy spectra were of interest from the point of 
view of investigation of quantum manifestations of classical stochasticity. 
Why should we address local characteristics of a spectrum? The fact is 
that the global characteristics such as the numbers of states N{E) or the 
smoothed density of levels p{E) are too crude. At the same time, such local 
characteristics as the function of the nearest-neighbour spacing distribution 
between levels (FNNS) are very sensitive to the properties of potential or to 
the shape of a billiard boundary. 
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6.1 The standard semiclassical spectrum un¬ 
folding method 


Already by definition statistical properties imply analysis of sufficiently large 
numbers of states, which necessarily lead us to the use of the semiclassical 
approach, which is well known to have as long and rich a history as quantum 
mechanics itself. 

The most fundamental characteristics of quantum spectra are the level 
density 

= ( 6 . 1 ) 

n 

and staircase state number function 

E 

n{E) = j dEp{E) = J2q{E- E^), (6.2) 


where E^ are stationary energy levels. 

In the semiclassical approximation it is common to consider separately so 
called smooth and oscillating components of the functions (6.1) and (6.2): 

p{E) = f,{E) + f,{E) 


n{E) = n{E) + n{E). 

The smooth component p{E) describes a gradual variation of the av¬ 
eraged levels density with energy growth, while the oscillating component 
describes the deviation of spectrum from the ” quasi-equidistant” one — the 
hypothetical spectrum with levels spacing equal to A{E) = l/p{E). At the 
present time there are very highly developed methods to obtain the semiclas¬ 
sical approximations (usually in the form of decomposition in powers of h) 
for p{E). Among them are the widely known Thomas-Fermi formulae, which 
are nothing but zero approximation for corresponding values: 

Ptf{E) = J dpdq6{E - H{p,q)) 

TiTFiE) = ( 2 ^ I dpdqQ{E - H{p, q)) 

where D is the dimensionality of the configuration space. As a rule in studies 
of quantum chaos we deal with Hamiltonians of the following type: 

2 

H{p,q) = ^ + U{q). 
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After integration on the momenta in (6.3), we obtain for the case D = 2 
(further we will limit ourselves to the two-dimensional systems only): 


Ptf{E) = ^ J dxdy6{E - U {x, y)) 
nTF{E) = ^ J dxdy{E - U{x, y))Q{E - U{x, y)) 


We can see that the Thomas-Fermi level density in two-dimensional dynami¬ 
cal systems is a non-decreasing function — constant for billiards and increas¬ 
ing for smooth potential systems. 


Quite often the Hamiltonian of the system is invariant to some discrete 
symmetry transformation, for example it has a parity of some kinc0 Cer¬ 
tain symmetry of the Hamiltonian allows us to reduce integration in ( |6.4[ ) 
over the whole {x, y) plane to integration over some small region, which may 
signihcantly simplify the task. It is even more important that, provided a 
certain symmetry of the Hamiltonian, the quantum states of the system di¬ 
vide on subsets — irreducible representations — according to their symmetry 
properties, for example even or odd states. States belonging to different ir¬ 
reducible representations are statistically independent. Therefore in order to 
research the statistical properties of such systems it is necessary to consider 
the states of different symmetry types separately. For that purpose it is nec¬ 
essary to know the semiclassical expressions for ’’partial” densities p^{E) and 
the number rR{E) of states, were A labels the symmetry type of the states. 
On that point there is an exhaustive general semiclassical theory |16], but 
in many concrete problems it is often easy to rediscover simple particular 
formulae. For example, in the simplest case of discrete symmetry — invari¬ 
ance with respect to reflection y —)■ —y — it is easy to obtain Thomas-Fermi 
expressions analogous to (6.3), for odd states separately. Indeed, the wave 
functions odd on y evidently satisfy the Schrodinger equation with potential, 
which differs from the original one only by the presence of an inhnitely high 
potential barrier situated in the y = t) plane. Effectively such impenetra¬ 
ble potential wall exactly halves the classically allowed phase space, and we 


^The presence of a continuous symmetry group in a two-dimensional dynamical system 
is of no interest, because it leads to the existence of an integral of motion additional to 
energy and therefore integrability of the problem excludes the possibility of dynamical 
chaos. 
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obtain: 


Ptf{E) = ^ f dxdyS(E -U{x,y)) 

y>0 

nTF(E)^^ f dxdp(E - U(x,p))e(E - U(x, p)) 


(6.5) 


y>0 


= ptf(E) - (^), = nTF(E) - n^/^(E). 

For zero — Thomas-Fermi — approximation we get a trivial result; 

p^j^^-(E) = p^^j^(E) = Ptf(E)/2 
'-(E) = n°r^(E) = nTF(E)/2. 


n. 


TF 


However such a symmetry decomposition approach is applicable in any order 
of ^-decomposition, where it already leads to non-trivial corrections. In the 
same way, but with much more difficulties, we can obtain the semiclassical 
approximations for partial level densities and staircase numbers with rota¬ 
tional symmetry — in Hamiltonian, invariant under rotations on ip = 2 f jn. 

Practically all two-dimensional models actively studied from the point 
of view of quantum chaology, namely the above mentioned potentials of 


quadrupole oscillation (2.3), Henon-Helis (3.3), lower umbilic catastrophe 


(2.4), coupled quartic oscillator (5.22) and also Barbanis potential 


UB(x,y) = 


x"^ + y"^ 


+ xy 


are described by rather simple expressions — polynomials of the third or 
fourth order. As a rule it allows us to analytically solve the equation U(x, y) = 
E with respect to, for example, variable y, and to obtain the explicit expres¬ 
sion for the level line y = f(x;E). It gives us the possibility to reduce 


the Thomas-Fermi integrals (6.5), and similarly the higher order semiclas¬ 


sical approximations, to one-dimensional integrals, which are evaluated if 
not analytically, then numerically in a much faster way than the original 
two-dimensional ones (Fig j6.1[ ). 

One of the most important applications of the above discussed semiclas¬ 
sical expressions is the so-called spectrum unfolding procedure using n{E): 


’^(En)- 

As a result the original spectrum is transformed to quasi-equidistant: 

p(e) = 1. 


( 6 . 6 ) 


Only for spectra transformed that way does it make sense to discuss the 
quantum signatures of chaos in the statistical properties of energy spectra. 
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Figure 6.1: Density p{E) (left) and staircase number n{E) (right) for h = 1 
in Thomas-Fermi approximation in the QO potential (2.3[). 


6.2 Nearest neighbor spacing distribution func¬ 
tion 

In the 1960s of last century Wigner m, Porter |1H] and Dyson lig built 
statistical theory of complex systems spectra on the basis of the following 
hypothesis: energy levels distribution in complex systems is equivalent to 
distribution of the eigenvalues of random matrices ensemble with certain 
symmetry. The ultimate result for FNNS obtained in that theory has the 
following form: 

P{S) ~ 

where a and b are slowly varying functions of energy and critical exponent 
a; dehning behavior of the distribution function in the limit S' —)■ 0, depends 
on symmetry properties of the random matrices: a = 1 for orthogonal and 
a = 2 for unitary ensemble of matrices. 

Predictions of the statistical theory of energy spectra (mainly for Gaus¬ 
sian orthogonal ensemble of random matrices) were carefully compared with 
all available data on nuclear spectra. No considerable contradictions between 
the theory and experiment were found. In particular, the random matrices 
ensembles perfectly reproduced such an important spectral characteristic as 
spectral rigidity, describing small fluctuations of energy levels around aver¬ 
aged values in a given interval. Measure of rigidity is the statistic A 3 of 
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Dyson |19] and Mehta [50] 


x-\-L 

A 3 (L; x) = y min j [n{e) — Ae — B]‘^de (6.7) 

L a,b J 

X 

which determines the least-square deviation of the staircase representing the 
cumulative density n{e) from the best straight line htting it in any interval 
[x,x + L\. The most perfectly rigid spectrum is the picket fence with all spac¬ 
ing equal (for instance, the one-dimensional harmonic oscillator spectrum), 
therefore maximally correlated, for which A 3 = 1 / 12 , whereas, at the oppo¬ 
site extreme the Poisson spectrum has a very large average value of rigidity 
(A3 = L/15), reflecting strong fluctuations around the mean level density. 

Analogous comparisons were made also for atomic spectra. Good agree¬ 
ment with theory was found for them too, although for much poorer statistics. 

A completely different approach to the problem of statistical properties 
of energy spectra was developed on the basis of non-linear theory of dynam¬ 
ical systems. As numerical simulations show ISHEa ESI EU, conhrmed by 
serious theoretical considerations [21 El Ea E6], the main universal property 
of systems which have regular type of dynamics in the classical limit, is the 
level clustering phenomenon, while for systems that are chaotic in the classi¬ 
cal limit level repulsion is observed. This statement sometimes is called the 
hypothesis of universal character of energy spectra fluctuations ISB. 

In the case of regular motion any classical trajectory lies on a surface 
topologically equivalent to a torus. Each such torus corresponds to a certain 
set of integrals of motions (or quantum numbers), with the total number 
equal to the number of degrees of freedom. Different eigenfunctions of inte- 
grable quantum systems correspond to different sets of quantum numbers, 
and therefore lie on different tori. Their eigenenergies are not correlated, 
which leads to Poissonian FNNS 

Pp[S) = e-^ ( 6 . 8 ) 


responsible for level clustering. 

In the transition to chaos some tori are destroyed, leading to formation 
of chaotic regions in the phase space. Classical trajectories in such regions 
diffuse between different (already destroyed) tori. Therefore different chaotic 
modes can be represented in the form of superposition of former regular 
modes, which leads to repulsion between the corresponding levels. In the 
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limit of a completely chaotic system, all modes of initially integrable sys¬ 
tem mix one with other, so that repulsion exists between any pair of levels, 
changing FNNS from Poissonian to the Wigner one: 

Pw{S) = (6.9) 


The situation became more complicated for generic Hamiltonian systems, 
where the phase space contains both regular and chaotic components. How 
is it reflected in the energy spectra statistical properties? 

Historically the hrst was the phenomenological approach to the solution 
of that problem, proposed by Brody m, who considered FNNS of the form: 


Ps^S) = (/? + 1)6 ^b 


F 



/ 3+1 


( 6 . 10 ) 


where (3 equals to the relative measure of classical phase space, occupied by 
the chaotic trajectories: /? = 0 for regular and /3 = 1 for chaotic systems. 

Berry and Robnik [58] and independently Bogomolny |59| using on semi- 
classical arguments, showed that FNNS in the case of mixed type dynam¬ 
ics represents the superposition of Poisson and Wigner distributions with 
weights that are determined by the relative phase space measure of regular 
and chaotic motion respectively: 


p _ d? 


e P^erfc(^{l—p)S 


p2e + |(1 — p) 25'^ e 4(1 


( 6 . 11 ) 

where p is the relative phase volume occupied by regular trajectories in the 
mixed system. The limit p —t 1 corresponds to regular system, and p —)• 0 
— to completely chaotic. 

Further development theory of FNNS acquired in the works of Narimanov 
and Podolskiy [60], who accounted for additional effect of repulsion for close 
levels due to dynamical tunneling (more details see Chapter]^: 


Pnp = P^F {^) e-P^erfc (^(1 - p),s) + 
+ (1 - p) (2pF (^) + 1(1 - pf_S) 

F{x) = 1 - 


( 6 . 12 ) 


e^—x 
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S 


Figure 6.2: The distribution functions of Poisson and Wigner (thin solid 
lines), Brody for = 0.5 (thick solid line), Berry-Robnik-Bogomolny for 
p = 0.5 (dashed line) and Narimanov-Podolskiy for p = 0.5 and v = 0.1 
(dotted line) 


The p parameter in the Narimanov-Podolskiy distribution function (6.12) 
has the same physical sense as in the Berry-Robnik-Bogomolny distribution 
(6.11), and additional parameter u is connected to the intensity of dynamical 
tunneling. Therefore the Berry-Robnik-Bogomolny distribution function rep¬ 
resents a particular case of the Narimanov-Podolskiy distribution at u = 9, 
which corresponds to the absence of the tunneling (Fig j6.2[ ). 

Among systems subject to extensive numerical analysis of the spectral 
properties, two-dimensional billiards possess the central place. The two- 
dimensional billiard represents a point particle freely moving on the plane in¬ 
side some region of arbitrary shape and elastically reflecting from the bound- 
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ary. Such systems attract active interest in research for the following reasons: 

1. the system has the lowest possible number of degrees of freedom, al¬ 
lowing chaotic motion in a conservative system; 

2. the simplicity of classical dynamics; 

3. stochasticity criteria for billiards are formulated in geometrical terms; 

4. homogeneity of the phase space; 

5. availability of efficient methods for solution of the Schrodinger equation 
for billiards; 

6. smooth component of level density is well known due to the Weyl for¬ 
mula; 

7. billiard dynamics reflect the real situation in many physical systems 
(quantum dots and antidots, Josephson’s junctions, nuclear billiards) 

8. at the present time there is efficient experimental realization both 
for classical and quantum dynamics in microwave and optical two- 
dimensional resonators. 

For billiards with certain boundary shape one of two limiting cases is 
realized: exact integrability or absolute chaos. So in the circular billiard 
(Fig j6.3^ ) angular momentum is the second (after energy) integral of mo¬ 
tion, and this system is integrable. A billiard of ’’stadium” type (Fig j6.3| 3) is 
one of the simplest stochastic systems. In Fig j6.3| the statistical characteris¬ 
tics of energy spectra (FNNS and dispersion) are presented for both systems. 
In complete agreement with the hypothesis of the universal character of spec¬ 
tral fluctuations, FNNS for the circular billiard is perfectly approximated by 
Poissonian distribution, and dispersion is a linear function of the length of 
the considered interval. In the non-integrable case the level repulsion effect 
is pronouncedly manifested, leading to Wigner distribution, and dispersion 
grows much more slowly because of the higher rigidity of the considered 
spectrum. 
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Figure 6.3: Statistical characteristics of energy spectra for circular (a) and 
’’stadium” (b) billiards. 
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6.2.1 Uniform potential case 


Let us note that an important advantage of billiards is the possibility of 


spectrum unfolding (6.6) with high precision and independence of relative 


measures of regular and chaotic components from the energy of particles. 
Both features equally share homogeneous potentials, in particular coupled 


quartic oscillators potential (5.22), where we consider for simplicity only 
Ai-type states — even on both coordinates and symmetric with respect to 
exchange x -ty y. 

As for quantum billiards, for the coupled quartic oscillators potential 


(5.22) there is a quite precise semiclassical expression for the smooth com¬ 


ponent of staircase state number function up to the third order in energy 




^,11 2 
^ f 2’2’^’~ 


a 


+ 6 


£(il 

23 


+ 


0+2 


r(|)y^ 


+ 


2Ei 


(6.13) 


which allows us to perform spectrum unfolding with sufficiently high accuracy 
in a wide energy range. However because of technical difficulties in numerical 


calculation of hypergeometric function in (6.13), in practice it is much more 


convenient to use the corresponding Thomas-Fermi approximation for (6.13) 
in the form: 




dy) 


1 + ^ sin^ 


(6.14) 


It is convenient to estimate the accuracy of the spectrum unfolding with a 
diagram of variable £„/n plotted for all levels in the spectral series under con¬ 


sideration. For a spectrum standardly unfolded by (6.6), this variable should 
be sufficiently close to unity in the whole energy interval. As can be seen 


from Fig, 6.4, the Thomas-Fermi approximation(6.14) gives slightly under 


stated values for fiA^i^E), but nevertheless allows us to achieve 1% accuracy 
of spectrum unfolding, and this accuracy improves with energy growth. 

In Fig j6.4| one can see gradual decrease of maximum amplitude of the 
spectral fluctuations with growth of chaoticity measure (see Fig 6.5), which 
already qualitatively points to gradual changes of the fluctuations character 
from Poissonian (Figj6.4^) to Wigner one (Figj6.4|i). 
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n 


n 


Figure 6.4: Accuracy estimation for the unfolding with (6.13) in the spectrum 
of coupled quartic oscillator potential (5.22) for a = 6(a), a = 7.5(b), a = 
11.5(c), a = 19.5(d). 
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The FNNS for coupled quartic potential ( |5.22 ) together with the best 
fits by Brody (6.10) and Berry-Robnik-Bogomolny (6.11) distributions are 
presented in Fig j6.6| in normal and on Fig j6.7| in logarithmic scale. The 
latter is especially convenient because in the logarithmic representation the 
Poissonian distribution turns into linear function: 


lnFp(^) = -^, 


and proximity of numerically obtained distribution to Poissonian one in the 
case of regular motion (Fig j6.7| ) becomes especially clear. 

According to the hypothesis of univeral character of spectral fluctuations 
the FNNS transforms from Poissonian at a = 6 to Wigner one at a = 19.5. 
In the intermediate cases the numerical distribution function is satisfacto¬ 
rily described by the Berry-Robnik-Bogomolny distribution as well as Brody 
distribution, and the latter gives even better accuracy regardless of its phe¬ 
nomenological nature. 

Distributions on Fig j6.6| and Fig j6.7| are obtained from analysis of spectral 
series consisted of 3000 levels each, but even these statistics appears too low 
to give preference to Brody distribution over the Berry-Robnik-Bogomolny 
distribution. In the case of too pure statistics which is common in smooth 
potentials it is better to use integral, or cumulative distribution: 


s 


hF(F) = 


P{S)dS 


0 


It is easy to obtain cumulative distributions corresponding to Brody and 
Berry-Robnik-Bogomolny distributions: 


kVERsiS) 


Wb{S) 
= 1 -1- A 

dS 


= 1 _ 


(6.15) 


Cumulative Poisson and Wigner distributions, as well as the usual ones, 
are limiting cases of the Brody distribution at /3 = 0 and /3 = 1 and of the 
Berry-Robnik-Bogomolny one at p = 1 and p = 0 respectively: 


WB{S) = l-e 


-s 




BRB 


{S) = l 


TT q 2 

e 4*^ 


(6.16) 


The numerically obtained cumulative FNNS for coupled quartic potential 


(5.22) together with the best hts by the Brody and Berry-Robnik-Bogomolny 
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Figure 6.6: FNNS for the coupled quartic oscillator potential (5.22) for a = 
6(a), a = 7.5(b), a = 11.5(c), a = 19.5(d). Points represent numerical data, 
solid lines — Poisson (6.8) and Wigner (6.9) distributions, dashed and dotted 
lines — the best hts by the Brody (6.10) and the Berry-Robnik-Bogomolny 
(6.11) distributions respectively. 
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(a) 


(b) 

0=7.5 





Figure 6.7: The same as on Figj6.6|but in logarithmic scale. 
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(6.15) distributions are presented in Fig j6.8[ It should be noted that param¬ 
eters of the best £t made independently by usual ( 6.10[6.11 ) and cumula¬ 
tive (6.15) distributions agree one with another very well. However cumu¬ 


lative distributions demonstrate better agreement with numerical data ex¬ 
actly for the Brody distribution rather than for the Berry-Robnik-Bogomolny 
distribution. The usual explanation of that fact states that the Berry- 
Robnik-Bogomolny distribution is a theoretically well-grounded universal 
high-energy asymptote for FNNS but it has limited applicability for insuffi¬ 
ciently high energies. Therefore level fluctuations in the spectral series avail¬ 
able for most potential systems are better described by phenomenological 


Brody distribution (6.10). 


The most striking distinction between the Brody and Berry-Robnik-Bogomolny 
distributions is displayed in the so-called T-representation: 

T(S) = ln{-ln[l-lT(^)]}. 

In the variables (T, InT) the Brody distribution with any parameter /3, 
as well as Poisson and Wigner distributions, become linear functions: 


Tb{S) = In 6+(/? + !) In 5 


(6.17) 


Tp(S) =lnS 
Tw(S) = In ^ + 2lnS, 

which is not the case for the Berry-Robnik-Bogomolny distribution; 


(6.18) 


Tbrb = In <j - In ( — 


e P^erfc(^^{l- p)S^ ^|. (6.19) 


The slope of the T-distribution on the (T, InT) determines the repulsion 


rate in the model (6.7): Wigner distribution corresponds to linear repul¬ 
sion and Poisson distribution — to zero, or absence of repulsion. Therefore 
to the Brody distribution we can assign so-called fractional levels repulsion 
with the exponent continuously changing from 0 to 1. At the same time 
the Berry-Robnik-Bogomolny distribution describes only integer repulsion 
with the exponent equal to 1. As can be seen from Fig, 6.9[ the fractional 
level repulsion effect actually takes place in smooth potential systems and is 
correctly described by the Brody distribution. 
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(a) 


(b) 



-f“=-1-1-'-1-'-1 -1-1-1-1-1-1 

0 1 2 3 0 1 2 3 


S S 


Figure 6.8: The cumulative FNNS for the coupled quartic oscillator potential 
(5.22) for a = 6(a), a = 7.5(b), a = 11.5(c), a = 19.5(d). Points repre¬ 
sent numerical data, solid lines — Poisson and Wigner (6.16) distributions, 
dashed and dotted lines — the best fits by the Brody and the Berry-Robnik- 
Bogomolny (6.15) distributions respectively. 
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(a) (b) 




(c) 



(d) 

0=19.5 



Figure 6.9: The cumulative FNNS in the T-representation for the coupled 
quartic oscillator potential (5.22) for a = 6(a), a = 7.5(b), a = 11.5(c), a = 
19.5(d). Points represent numerical data, solid lines — Poisson and Wigner 
(6.18) distributions, dashed and dotted lines — the best hts by the Brody 
(6.17) and the Berry-Robnik-Bogomolny (6.19) distributions respectively. 
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6.2.2 Generic potential case 


Let us pass now from the special case of uniform potentials to the more 
generic situation. Let us note at the beginning that, even in the case of non- 
uniform potentials in the presence of R — C or R — C — R transitions, it is 
possible to build rather representative series of energy levels corresponding to 
a dehnite type of classical motion. As an example let us consider the trans¬ 
formation of statistical properties of the energy spectrum of the quadrupole 
oscillations Hamiltonian (2.2) in the R — C — R transition for one-well case 


{W < 16). At the hxed topology of potential surface (W = const) the unique 
free parameter of the Hamiltonian is the scaled Planck’s constant h. In the 
study of the concrete energetic interval (i?i, C or R 2 ), corresponding to a 
dehnite type of classical motion, the choice of h is dictated by the possibility 
of attainment of the necessary statistical assurance (the required number of 
levels in the investigated energy interval) with conservation of precision of 
spectrum calculation (restrictions to possibility of diagonalization of matrices 
of large dimension). The numerical results are presented in Fig.6.10 Both 


the FNNS p{x) and the average value of A 3 well correspond to the predic¬ 
tions of GOE for the chaotic {C) region. The logarithmic scale for p{x) is 
suitable to trace this correspondence at large x. For regular regions (i?i and 
R 2 ) the distribution function, in the same scale, according to the hypothesis 
of the universal character of huctuations of energy spectra, must be repre¬ 
sented by a straight line (the logarithm of Poisson’s distribution). The results 
demonstrate the agreement with this hypothesis, though small-sized devia¬ 
tions are observed for small distances between levels. Such a tendency to the 
rise of some repulsion in the regular region, apparently, is connected with 
a small admixture of chaotic components. At the construction of statistical 
characteristics, a purity of sequence is provided by using only those levels, 
which are relative to dehnite irreducible representation of group (the 
levels of E-type were used for results represented in Fig jG.lOt the statistical 
characteristics of levels of Ai and A 2 -types have similar form). 

Spectral ranges of multi-well potentials realizing the mixed state open new 
possibilities for investigation of the intermediate statistics. At those energies 
chaotic and regular components are separated not only in Phase space as in 
the usual case of mixed type motion), but already in the conhguration space. 
A priori the FNNS form in the mixed state is not necessarily reduced to 
dehnitely weighted superposition of Poisson and Wigner distribution. In that 
case we deal not with statistics of admixture of two level series with different 
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Figure 6.10: Correlation between the character of classical motion and statis¬ 
tical properties of energy spectra in R — C — R transition for the Hamiltonian 

On the left — Poincare sections, in the mid- 
the variance T?. From the 


of quadrupole oscillations (2.2). 
die 


FNNS distribution p(a;), on the right 
bottom up: the hrst regular range i?i, the chaotic range (7, the second regular 
range i ?2 
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Figure 6.11: FNNS for the potential (2.4) for E < Ecr{a) and E > 
Ecr{h). Points represent numerical data, solid lines — Poisson (6.8) and 
Wigner ( |6.9 ) distributions, dashed line — the best £t by the Berry-Robnik- 
Bogomolny (6.11) distribution. 


FNNS, but with statistics of spectrum series, where each level cannot be 
attributed to definite (Poisson or Wigner) statistics. Statistical properties of 
such systems have not been under study up to now although they are the 
systems corresponding to the common case situation. 

Let us consider the simplest potential realizing the mixed state — two- 
well potential of lower umbilic catastrophe (2.4). In order to describe 
the statistical properties of the corresponding energy spectrum, let us try to 
use the Berry-Robnik-Bogomolny distribution (6.11), where the p parameter 
equals the relative phase volume occupied by regular trajectories. In the 
case of potential p ~ 1 for E < E^r and p ~ 1/3 for E > E^r, which 
qualitatively agree with numerical data for energy levels statistics Fig j6.11| 

The Berry-Robnik-Bogomolny distribution function by its construction 
describes the energy level fluctuations in the absence of interaction between 
the regular and chaotic components. Therefore the best agreement between 
the Berry-Robnik-Bogomolny spectral fluctuations theory and experimental 
data we should expect exactly in the mixed state spectra, where interaction 
of chaotic and regulat states each with other is additionally suppressed by 
the energy barrier separating the corresponding local minima. In practice 
however such agreement is never observed for several reasons. 

The hrst reason is not connected immediately with the mixed state prop¬ 
erties and lies in technical difficulties of calculation of long spectral series in 
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potentials with several local minima. In practice one usually should be happy 
with a series of several hundred levels at most, which is dehnitely insufficient 
to obtain accurate distribution functions for spectral fluctuations. Another 
difficulty of a technical nature is connected with the fact that the multiplic¬ 
ity of generic smooth potentials does not allow any systematic semiclassical 
analysis, because in particular there are no sufficiently precise semiclassical 
approximations for the smoothed component of staircase states number func¬ 


tion, required for the standard unfolding procedure (6.6). Therefore in order 


to unfold the spectra we should be happy with the zero semiclassical approx¬ 
imation, i.e. usual Thomas-Fermi approximation, which for even states in 
the potential takes the following form: 


nR,c{E) 


n{E) = nR{E) + nc{E) 

^/2(l+^/E) 

/ 


y/2{l-CE) 


dx 

/^±2 


3 




( 6 . 20 ) 


where nR{E) and nc{E) represent partial contributions from regular and 
chaotic minima respectively. Such a simple approximation nevertheless al¬ 
lows us to perform the spectrum unfolding with maximum errors of order 1% 
in the whole energy interval of interest (see Fig j6.12[ ). 

Therefore in potentials with multiple local minima, the smoothed semi¬ 
classical state number function can be calculated as the sum of partial con¬ 
tributions from each local minimum. In the case of mixed state it gives 
a convenient way to estimate the relative density of states of regular and 
chaotic type. 

The second reason for disagreement between Berry-Robnik-Bogomolny 
distribution and actual numerical data is common for all potential systems 


(excluding uniform potentials like the coupled quartic oscillator one (5.22)), 


and it consists of the fact that in smooth potential systems the classical 
chaoticity measure essentially changes with energy; therefore upper and lower 
levels in sufficiently long spectral series will correspond to different relations 
of regular and chaotic components, while all known theoretical distribution 
functions assume the relation to be constant. It forces us to consider spectral 
series that are sufficiently narrow in energy, and this makes levels statistics 
even poorer. 

Finally, the specific feature of potential, which it shares with D^, 
Barbanis and Henon-Heilis potential, consists of the fact that in all the 
above mentioned potentials the discrete energy spectrum, strictly speaking. 











104 


V.P.Berezovoj, Yu.L.Bolotin, V.A.Cherkaskiy, G.I.Ivashkevich 


(a) 


0.50<E<0.56 


(b) 


0.60<E<0.66 








(c) 


0.93<E<0.96 


(d) 


0.96<E<0.99 



n 


Figure 6.12: Accuracy estimation for the unfolding with (6.20) in the spec¬ 
trum of lower umbilic catastrophe potential (2.4) for 0.50 < E < 0.56(a), 
0.60 <E < 0.66(b), 0.93 <E < 0.96(c), 0.93 <E < 0.99(d). 
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is absent. It is due to the fact that all those models allow inhnite motion 
for all energies. Nevertheless many researchers discuss spectral fluctuations 
in such systems, implying the spectrum of quasistationary states, localized 
in the corresponding potential well and extremely slowly decaying due to 
tunneling into the continuous spectrum. Any numerical calculation of such 
states practically contains implicit reformulation of the original model in 
order to make it more correct from the point of view of quantum mechanics. 
The quasi-stationarity of such states becomes especially manifest near the 
saddle energies — the most interesting energy region for investigation of the 
mixed state. 

The FNNS for potential obtained for narrow energy intervals are 
presented in Fig, 6.13| in normal and in Fig, 6.14 in logarithmic scale, and 
the corresponding cumulative distributions are presented in Fig j6.15| inlF- 
representation and on Fig j6.16 in T-representation. Regardless of the rather 
poor statistics (only 100 levels in each spectral series), we should note quite 
good qualitative agreement between the htting parameters of the Brody 
and the Berry-Robnik-Bogomolny distributions on the one hand and rela¬ 
tive phase volumes of regular and chaotic regions in classical phase space 
(see the Poincare surfaces of section on Fig j6.17[ ). 

As we can see, the Brody distribution gives better agreement than the 
Berry-Robnik-Bogomolny distribution with numerical FNNS distributions in 
the potential. It is manifested both in the accuracy for each energy 
interval separately and in the variation of (3 parameter with energy growth; 
in contrast to the p parameter in the htting by Berry-Robnik-Bogomolny 
distribution, it grows monotonously with growth of chaoticity measure of 
classical motion and coincides in order of magnitude with the relative volume 
of the classical chaotic component. 

The unsatisfactory applicability of Berry-Robnik-Bogomolny the distri¬ 
bution for description of spectral huctuations in the mixed state is explained 
hrst of all by insufficient quasiclassicality of the considered states. However 
a more important reason for such disagreement of theory and experiment 
appears to be connected with the specihc features of energy spectra in the 
mixed state. It lies in the fact that that many of the energy levels in the po¬ 
tentials with multiple local minima actually appear to be very close to those 
energy levels that correspond to each local minimum separately — neglect¬ 
ing the small inhuence of other minima. Those levels are in their turn very 
well described by spectral series obtained in the harmonic oscillator approx¬ 
imation, especially for not too high energies. In the case of the mixed state 
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Figure 6.13: FNNS for lower umbilic catastrophe potential (2.4) 

for 0.50 < E < 0.56(a), 0.60 < E < 0.66(b), 0.93 < E < 0.96(c), 
0.93 < E < 0.99(d). Points represent numerical data, solid lines — Poisson 
(6.8) and Wigner (6.9) distributions, dashed and dotted lines — the best hts 
by the Brody (6.10) and the Berry-Robnik-Bogomolny (6.11) distributions 
respectively. 





















Regular and Chaotic Classical and Quantum Dynamics 


107 



(C) 



(d) 



Figure 6.14: The same as on Figj6.13|but in logarithmic scale. 
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(a) (b) 




Figure 6.15: Cumulative FNNS for lower umbilic catastrophe potential 
Q for 0.50 < E < 0.56(a), 0.60 < E < 0.66(b), 0.93 < E < 0.96(c), 
0.93 < E < 0.99(d). Points represent numerical data, solid lines — Poisson 
and Wigner ( |6.16| ) distributions, dashed and dotted lines — the best hts by 
by Brody and Berry-Robnik-Bogomolny (6.15) distributions respectively. 
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(a) (b) 




(c) (d) 




Figure 6.16: Cumulative FNNS in the T-representation for lower umbilic 
catastrophe potential (2.4) for 0.50 < E < 0.56(a), 0.60 < E < 0.66(b), 
0.93 < E < 0.96(c), 0.93 < E < 0.99(d). Points represent numerical data, 
solid lines — Poisson and Wigner (6.18) distributions, dashed and dotted 
lines — the best hts by the Brody (6.17) and the Berry-Robnik-Bogomolny 
(6.19) distribntions respectively. 
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Figure 6.17: Poincare surfaces of section for lower umbilic catastrophe 
potential (2.4) for E = 0.5(a), E = 0.6(b), E = 0.9(c), E < 1(d). Solid lines 
limit classically allowed region of the phase space. 
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especially close to the harmonic oscillator levels appear to be those states 
that correspond to minima with regular motion in the classical limit. 

In potential harmonic approximation is given by a mixture of two 
spectral series of two-dimensional harmonic oscillator levels with eigenfre- 
quencies 

Ux — 2 


2^/1 


where the eigenfrequencies are different for regular and chaotic local 


minima. 


Fig,6.18 demonstrates unexpectedly good agreement (in the limits 


of 1% of mean level spacing for not too high energies) between the harmonic 
approximation and actual energy levels in the potential, and this agree¬ 
ment is better for the regular minimum. 

Such proximity of the considered spectrum to the harmonic one automati¬ 
cally implies also the mutual proximity of corresponding FNNS distributions. 
As is well known, FNNS distribution for the harmonic oscillator has a patho¬ 
logical nature — it does not have any universal form, and, therefore, is not 
described by any of the known distribution functions. Fig j6.19j|6.23| present 
the results of analysis of spectral fluctuations in the harmonic oscillator ap¬ 
proximation in lower umbilic catastrophe potential ( 2.4[ ) for chaotic and 
regular minima separately, and also for a spectrum obtained by mixing of 
spectral series of both minima. 

The FNNS for harmonic oscillator spectrum has multiple maxima, de¬ 
termined by its eigenfrequencies. Therefore the FNNS distribution for the 
harmonic approximation spectrum in the mixed state has many maxima, 
corresponding to all eigenfrequencies of the approximations. That depen¬ 
dence of the FNNS distribution form on the eigenfrequencies of the system 
is the main reason for the pathological character of level fluctuations in the 
harmonic oscillator spectrum. In the case of the mixed state, energy levels, 
close to the harmonic ones, will give essential (and pathological) contribution 
to FNNS distribution. The theoretical distribution function, accounting for 
such a contribution, should be constructed as mixed statistics from compo¬ 
nents of three types: not only chaotic or regular, but also harmonic ones, 
carrying pathological spectral fluctuations. The construction of such a dis¬ 
tribution function is a problem for the future. 









112 


V.P.Berezovoj, Yu.L.Bolotin, V.A.Cherkaskiy, G.I.Ivashkevich 



n 


(b) (c) 




Figure 6.18: Energy levels of lower umbilic catastrophe potential (2.4): 
(a) numerically exact spectrum (solid line), harmonic oscillator (dashed line) 
and the semiclassical (dotted line) approximations; (b) absolute value of 
difference between the numerically exact energy levels and harmonic oscil¬ 
lator approximation; (c) absolute value of deviance of numerically exact en¬ 
ergy levels from the harmonic oscillator approximation referred to mean level 
spacing. 
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Figure 6.19: Accuracy estimation for the spectrum unfolding of harmonic 
approximation in lower umbilic catastrophe potential (2.4): (a) for chaotic 
minimum {x < 0 ); (b) for regular minimum [x > 0 );(c)for both minima 
together. 



Figure 6.20: FNNS for harmonic approximation in lower umbilic catastrophe 
ZI 5 potential ( 2.4[ ): (a) for chaotic minimum (x < 0); (b) for regular minimum 
{x > 0);(c)for both minima together. Points represent numerical data, solid 
lines — Poisson ( 6 . 8 ) and Wigner (6.9) distributions. 



Figure 6.21: The same as Fig,6.20 but in logarithmic scale. 
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Figure 6.22: Cumulative FNNS for harmonic approximation in lower umbilic 
catastrophe potential ( 2.4[ ): (a) for chaotic minimum (x < 0); (b) for 
regular minimum (x > 0);(c)for both minima together. Points represent 
numerical data, solid lines 

(a) 


Poisson and Wigner (6.16) distributions. 

(b) - (c) 



Figure 6.23: Cumulative FNNS for harmonic approximation in lower um- 
bilic catastrophe potential (2.4) in the T-representation: (a) for chaotic 
minimum (x < 0); (b) for regular minimum (x > 0);(c)for both minima to¬ 
gether. Points represent numerical data, solid lines — Poisson and Wigner 
(6.18) distributions. 
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6.3 Other statistical characteristics of spec¬ 
tra — l//-noise 


FNNS is far from being the unique tool for investigation of quantum manifes¬ 
tations of classical stochasticity in the structure of quantum energy spectra 
in the systems whose classical analogues demonstrate chaos. One of the 
alternative approaches [61] to quantum chaos description is based on tradi¬ 
tional time series analysis methods. Let us consider the energy spectrum as 
a discrete signal, and a sequence of N energy levels as a time series. Let 
us characterize the spectral fluctuations with the ^^-statistics, dehned as the 
following: 

n n 

<5n = - (s)) = 

i=\ i=l 


where index n runs from 1 to — 1. Quantities Wi give the fluctuation of 
the f-th spacing between neighboring levels with respect to its mean value 
(s) = 1. We assume that the standard procedure of spectrum unfolding (6.6) 
was made in advance, i.e. the original spectrum was mapped into another 
spectrum (L’j —>■ Si) with mean level spacing equal to unity, then Sj = Sj+i—£j. 


i = 1,... ,iV. 

Our goal is to study transformation of the above introduced characteris¬ 
tics 5n in transition from regular systems to chaotic ones. One uf the useful 
methods is to calculate the power spectrum S{k) of the discrete and hnite 
series 5^ according to the following: 


S{k) = |4|2 


where 5k is the Fourier transform of (5„, 

1 


5k. = 


^/N 




and N is the number of terms in the series. 

As an example of a chaotic system authors [61] have chosen the atomic 
nucleus at high excitation energies, where level density is very high. In or¬ 
der to obtain the energy spectrum the authors performed calculations in the 
shell model with realistic interaction reproducing the experimental data well. 
Hamiltonian matrices for differen values of angular momentum were diago¬ 
nalized and the unfolding operation was accurately done for the obtained 
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spectra. After that the authors selected sets from 256 consecutive levels each 
with the same J^, T from high level density regions. In order to characterize 
the statistical properties of the signal Sn the authors calculated the ensemble 
average of its power spectrum. The averaging was performed in order to de¬ 
crease the statistical fluctuations (errors) and to clarify the general tendency. 
Average {S{k)) was calculated over the ensemble of 25 sets. 

Fig,6.24 represents the results for a typical stable sd shell of ^^Mg nucleus 


(the spectrum was obtained by diagonalization of a 2000 x 2000 martrix) and 
for a very exotic ^^Na nucleus — sd proton shell and pf neutron one (5000 x 
5000 matrix diagonalization). It can be clearly seen that the power spectrum 
for Sn is close to power law. It is possible to assume simple functional form 


(^(^)) - 


The least squares £t gives a = 1.11 ± 0.03 for ^'^Na a = 1.06 ± 0.05 
for ^“^Mg. A natural question arises: does there exist an universal relation 
between peculiarities of quantum spectrum, connected with the character of 
classical motion, and the power spectrum of Sn fluctuations? 

Probably the simplest way to answer that question is to compare Sn and 
{S{k)) for level sequences corresponding to Poissonian statistics with the 
results for different random matrices ensembles. Fig, 6.25 presents the Sn 
signal for the Poissonian and GOE spectra of dimensionality 1000. It is clear 
that the signals are absolutely different. We compare them with discrete time 
series x{t) corresponding to power spectra 1/k and l/fc". The proximity of 
the a = 2 time series to the Poissonian spectrum and of a = 1 — to to GOE 
one is evident. 

In order to calculate {S{k)) there were built 30 different matrices of di¬ 
mensionality 1000 X 1000 for random matrix ensembles of different types 
mentioned above. Fig, 6.26 presents the results of the calculations in double 
logarithmic scale. In all the considered cases the general tendency is linear 
with the exception of the very high energy region, where some deviations 
are observed which are probably due to the hniteness of dimensionality of 
the considered matrices. Ignoring frequencies with k > 2.2, the htting gives 
a = 1.99 for Poissonian spectrum with uncertainty about 2%. 

In contrast to the last example, the GOE of high dimension matrices 
is usually considered as a paradigm for the chaotic quantum spectrum. It 
demonstrates level-level correlations on all scales. The same refers also to 





loR<S(k)> 
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Figure 6.24: Average power spectrum of the 6n function for ^^Mg and ^"^Na, 
using 25 sets of 256 levels from the high level density region. The plots are 
displaced to avoid overlapping. 
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Figure 6.25: Comparison of the 5n function for Poisson (dashed line) and 
GOE spectra (solid line), with a standard time series x{t) with l//c“ power 
spectrum, for a = 2 (dashed line) and a = 1 (solid line). 
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Figure 6.26: Power spectrum of the 6n function for GDE (Poisson) energy 
levels, compared to GOE, GUE and GSE. The plots are displaced to avoid 
overlapping. 


unitary (GUE) and simplectic (GSE) ones, but with increasing level repulsion 
rate. 

Power spectra of 6n for all the three ensembles are presented on Figj6.26 


For the exponents we obtained the following values: acoE = 1-08, acuE = 
1.02, acsE = 1.00. All three ensembles agree with the same power law a ~ 1. 
It is evident that the power spectrum {S{k)) behaves as l/fc" both for regular 
and chaotic spectra, but level correlations decrease from limiting a = 2 for 
a regular uncorrelated spectrum to the minimal value a = 1 for quantum 
chaotic systems. 

The above result for the 5n statistics power spectrum gives rise to the 
following hypothesis: energy spectra of chaotic systems are characterized by 
1// noise. Let us recall that 1// (flicker noise) was discovered in the 1920s: 
slow flicks (from which the name comes) of emission capacity of cathode 
lamps (as a result of current fluctuations) superposed on fast thermal (frac- 
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tion) current fluctuations. For white noise S{uj) = const in a wide frequency 
range, while for the flicker noise *S'(a;) (x 1/uj. 

Flicker noise was discovered in many physical and biological systems. A 
common property of all such systems is the fact that 1// noise is a satellite of 
stationary irreversible processes: its contribution disappears in the absence 
of irreversible flows and the system reaches the thermodynamic equilibrium. 
In reality 1// noise is always weak in the sense that it carries only a small 
part of average fluctuation squared: the main contribution is due to white 
noise. But at low frequencies 1// noise may exceed white noise by several 
orders of magnitude and becomes the main harmful factor for devices that 
work on such frequencies. An important feature of flicker noise is the increase 
of its spectral density with frequency decrease, which takes place up to min¬ 
imal measurable frequencies of about 10“^ Hertz without any tendency to 
saturation. Thus the main question is: should one seek for a universal flicker 
nose mechanism common for all media and systems, or it is better to assume 
the existence of many flicker noises of different characters. This question is 
unanswered up to now. 

The hypothesis proposed in [61] has several attractive features. The con¬ 
sidered properties characterize the chaotic spectrum immediately without 
reference to properties of other systems (such as random matrices ensem¬ 
bles). It is universal for all types of chaotic system and does not depend 
on such properties as time reversal invariance and whether spins are integer 
or half-integer. Also the 1// characteristics unites a rich diversity of quan¬ 
tum chaotic systems from many different fields of science — flicker noise is 
everywhere! Therefore the energy spectra of chaotic systems demonstrate 
the same type of fluctuations as many other systems. However there is no 
indication that 1// noise in spectral fluctuations of quantum systems implies 
1 // noise in the corresponding classical analogues. 


Chapter 7 

Signatures of Quantum Chaos 
in Wave Function Structure 

7.1 QMCS in wave functions structure 

As we have seen in the previous section, the statistical properties of spectra 
have been rigidly correlated with the type of classical motion. It is natural to 
try to discover the analogous correlations in the structure of wave functions, 
i.e. to assume that the form of the wave function for a semiclassical quantum 
state depends on whether it is associated with classical regular or chaotic 
motion. Furthermore it should be pointed out that in the analysis of quantum 
manifestations of classical stochasticity at the level of energy spectra, the 
principal role was given to statistical characteristics, i.e., quantum chaos was 
treated as a property of a group of states. The choice of a stationary wave 
function of the quantum system, which is chaotic in the classical limit, as 
the basic object of investigation relates the phenomenon of quantum chaos 
to an individual state. Joint investigation of both the possibilities is not a 
contradiction but reflects the universality of the considered phenomenon. 

In contrast to spectrum, the form of wave functions depends on the ba¬ 
sis on which they are determined. For QMCS studies the following three 
representations are most often used: 

1. The so-called FTo-i'epi'esentation is the representation of eigenfunctions 
{(Pn\ of the integrable part Hq of total Hamiltonian H = Hq + V. The 
main objects of investigation in this case are the coefficients of expan¬ 
sion Cmn of stationary functions 'ipm on basis iJo-representation 
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is natural in the realization of numerical calculations, as the diagonal- 
ization of Hamiltonian H is realized more often in this representation. 

2 . Coordinate representation \h(x) is the most convenient from the point 
of view of visual clarity as it allows direct comparison with classical 
motion along the trajectories in coordinate space. Indeed, according to 
Schnirelman theorem [02], the mean probability density |'h(x)p in the 
semiclassical limit h —)■ 1 coincides with projection of micro-canonical 
distribution onto the coordinate space: 

(|\l'(g)P) —)■ po{q) with h —)■ 1 

j(Bp5{E-H{p,q)) 

/ d^p'd^q'6{E - H{p', q')) 

Actually this theorem transmits the correspondence principle on wave 
functions. What can we expect from the wave function structure bas¬ 
ing on the theorem? It is concentration of the probability density in 
the regions of coordinate space covered by (quasi-)periodic trajectories 
for regular wave functions, in contrast to almost uniformly smeared 
probability density for chaotic wave functions. 

3. Representation with the help of Wigner functions [HSl El] has a set 
of properties in common with the classical function of distribution in 
phase space. 

As early as 1977 Berry [HS] assumed that the form of the wave function 
'ip for a semiclassical regular quantum state (associated with classical motion 
on a Wdimensional torus in the 2 Wdimensional phase space) was very dif¬ 
ferent from the form of pj for an irregular state (associated with stochastic 
classical motion on all or part of the {2N — l)-dimensional energy surface in 
phase space). For regular wave functions, the average probability density in 
the conhguration space was seen to be determined by the projection of the 
corresponding quantized invariant torus onto the configuration space, which 
implies global order. The local structure is implied by the fact that the wave 
function is locally a superposition of a finite number of plane waves with the 
same wave number as determined by classical momentum. In the opposite 
case for chaotic wave functions the averaged (over small intervals of energy 
and coordinates) in the semiclassical limit h —)■ 0 coincides with the pro¬ 
jection of the classical micro-canonical distribution onto the coordinate space. 
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Its local structure is spanned by the superposition of infinitely many plane 
waves with random phases and equal wave numbers. The random phases 
might be justihed by classical ergodicity, and this assumption immediately 
predicts locally the Gaussian randomness for the probability of amplitude 
distribution. Such a structure of wave function is in good agreement with 
the picture of classical phase space: the classical trajectory homogeneously 
hlls the isoenergetic surface in the case of chaotic motion. By contrast, from 
the consideration of regular quantum state as an analogue of classical mo¬ 
tion on a torus, a conclusion should be done about the singularity (in the 
limit h —)■ 0) of the wave function near caustics (boundaries of the region of 
classically allowed motion in a coordinate space). 

Berry’s hypothesis was subjected to the most complete test for billiards of 
different types and, in particular, for a billiard stadium type |SH]- The ampli¬ 
tude of a typical wave function of an integrable circular billiard is negligible 
in the classically prohibited region (from conservation of angular momen¬ 
tum in a circular billiard it should be clear that the arbitrary trajectory 
is enclosed between the external and a certain internal circle, the radius of 
which is determined by the angular momentum magnitude) , whereas near 
the caustics it is maximal. Distribution of iV’nP for the case of a stadium 
billiard, (classical dynamics is stochastic) is strongly distinguished from the 
integrable case. However, this distribution is not so homogeneous, as we 
could expect starting from the ergodicity of classical motion |S7|. 

The above considered R — C — R transition represents an attractive pos¬ 
sibility to check Berry’s hypothesis for potential systems with a localized in¬ 
stability region. Let us start from the iLo-representation, or more exactly the 
representation of linear combinations of wave functions of a two-dimensional 
harmonic oscillator with equal frequencies 

‘Pi = J2<^’-NLi\NLj), ( 7 . 1 ) 

N,L 


where 


\NL3) = :^(|iV,L) +j|iV,-L)), j = ±1; 
AT = 0,1, 2,...; L = AT, iV - 2,..., 1 or 0; Pyj = 




{NLj\N'L'f) 
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Figure 7.1: Distribution of the coefficients on the number i = {A^, L,j} 
of the basis states. The numbers k of the corresponding levels are shown over 
the figures. 


If we introduce the notion of distributivity of wave function in basis, then 
the criterion of stochasticity formulated by Nordholm and Rise |6H] states 
that the average degree of distributivity of wave functions arises with the 
degree of stochastization in the system. It is clear that this criterion is a di¬ 
rect analog of Berry’s hypothesis for Ffo-representation, if one interprets the 
number of the basis state i^ij as a discrete coordinate. Fig ]7.1| qualitatively 
conhrms this criterion. It can be seen from this hgure that the states that 
correspond to regular motion (regions Ri and R 2 ) are distributed in a rela¬ 
tively small number of basis states. At the same time, states corresponding 
to chaotic motion (region C) are distributed in a considerably larger number 
of basis states. In the latter case, the contributions from a large number of 
basis states in expansion (7.1) interfere; this results in a complicated spatial 
structure of the wave function %lj[x,y). 

Correlations between the structure of wave function and the type of clas¬ 
sical motion are also demonstrated in hg j7.2[ where the probability density 
\'^k{,x,y)\^ for states with numbers 111, 160 and 210 is represented [9]. The 
squared module of the wave functions reproduces rather well the transition 
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from functions with a clear internal structure (region RQ to an irregular dis¬ 
tribution (region C) and the restoration of structure in the second regular 
region {R 2 )- For the chosen technique, in which transition is traced for the 
wave function with a hxed number ( scaled Planck constant), a change in the 
wave function is associated only with R — C — R transition. 

The specihc character of stationary wave functions corresponding to a 
dehnite type of classical motion manifests also in the level lines structure, 
in particular in the structure of manifolds T(q) = 0. Depending on the 
conhguration space dimensionality, those manifolds are called nodal points 
{N = 1), lines {N = 2), or surfaces {N >2). 

Already in 1979, Stratt, Handy and Miller [69] had proposed a quantum 
chaos criterion directly connected with nodal geometry: for regular states 

the system of nodal lines T(q) = 0 represents a grid of quasi-orthogonal 
curves (or is very close to it),while for chaotic states there is no such repre¬ 
sentation. 

For two-dimensional systems it is more convenient to describe not the 
nodal lines themselves, but the regions of constant sign of the wave function 
— the so-called nodal domains, and then the borders between the nodal 
domains coincide with the nodal lines. As can be easily seen from hg j7.3[ the 
stationary wave functions nodal structure in integrable systems may have 
properties usually assigned to chaotic ones — it is an irregular picture of 
nodal domains and avoids intersections of nodal lines. 

Fig j7.4| conhrms that the structure of the lattice of nodal curves of QO 
wave functions undergoes a change in the R — C — R transition. The spatial 
structure of nodal curves for states from regions Ri and i ?2 of regular classical 
motion is considerably simpler than this structure for states from the chaotic 
region C. 

The considered possibility corresponds to the traditional approach in 
search of QMCS in wave function structure: it is the investigation of changes 
in the wave functions structure with energy increasing. Such an approach 
adresses the difficulties connected with the necessity to separate the QMCS 
from the trivial modihcations in the wave functions structure due to quantum 
numbers changes. 

Let us now turn to the multi-well case in the potentials, for example, 
QO with hF > 16 and D^. The stationary wave functions corresponding 
to the energies at which the classical mixed state is observed, allow us to 
see the QMCS in comparison of different parts of one and the same wave 
function [Toiin]. As in the mixed state, different types of classical motion 
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k=111 k=160 



k=213 










Figure 7.2: Isolines of probability density The step between the 

isolines is constant. The numbers k of the corresponding levels are shown 
over the hgures. 
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Figure 7.3: The nodal domain structure for the rectangular billiard stationary 
wave functions ^(a;, y) = sin Stto; cos Stt?/ + £ sin Stto; cos Stt?/ for e = 0.01 (a), 
e = 0.3 (b), e = 0.6 (c) and e = 0.99 (d). 







Figure 7.4: The nodal curves for the wave functions \hfc(a;,?/). The numbers 
k of the corresponding levels are shown over the hgures. 













Regular and Chaotic Classical and Quantum Dynamics 


129 


are realized in different local minima; it is manifested in the wave functions 
structure. The parts of a wave function corresponding to different local min¬ 
ima show the above described features, characteristic for regular or chaotic 
motion respectively. Comparing the structure of the eigenfunction in central 
and peripheral minima of the QO potential (Fig j7.5^ ,c), or in left and right 
minima of the potential (Fig‘i ]7.5| 3,d), it is evident that the nodal struc¬ 
tures of the regular part and the chaotic part of the eigenfunction are clearly 


different. 


i) within the classically allowed region the nodal domains of the regular part 

of the wave function form a well recognizable checkerboard-like pattern 
[ST], while nothing similar can be observed for the chaotic part; 

ii) the nodal lines of the regular part exhibit crossings or very tiny quasi¬ 

crossings; in the chaotic part the nodal lines quasi-crossings have sig- 
nihcantly larger avoidance ranges; 

iii) while crossing the classical turning line U{x,y) = Em the nodal lines 

structure of the regular part immediately switches to the straight nodal 
lines, going to inhnity, which makes the turning point line itself easily 
locatable in the nodal domains structure; in the chaotic part, an inter¬ 
mediate region exists around the turning line, where some of the nodal 
lines pinch-off, making transition to the classically forbidden region 
more gradual and not so manifest in the nodal structure. 

Let us note an additional advantage in nodal structure investigations 
for QMCS search in the mixed state case. In generic multi-well potentials, 
stationary wave functions as a rule are strongly localized in one or a few 
local minima. The localization rate is conveniently characterized by the 
dimensionless variable (for example in the case of two-well potential DQ 


pR _ pC 
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^ n. ' ^ n. 


n ' n 
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Figure 7.6: Localization of the wave functions in potential. 


Evidently p = 1 for regular and r] = —1 for chaotic states. Proximity of 
Tj to unity is characterized by the variable 

U = -sign{r]n) logio(l - \Vn\) 

Therefore = 1 corresponds to a wave function for 90% localized in the 
regular minimum, and — 1 — for 90% in the chaotic one. In the mixed 

state case, the amplitude of overwhelming majority of wave functions in 


different local minima differs by several orders of magnitude (£g, 7.6), which 
makes difficulties for analysis of the probability density distribution, but 
those difficulties are however inessential for the analysis of wave functions 
nodal structure. 


7.2 Evolution of shell structure in the process 
oi R — C — R transition 

Introduction of the stochasticity concept in nuclear theory made it possible 
to take a fresh look at the old paradox: how one could reconcile the liquid 
drop — short mean free path — model of the nucleus with the independent 
particles — gas-like shell model. To solve this paradox it is sufficient to 
assume the following: 

1. When the nucleonic motion inside the nucleus is integrable, we expect 
to see strong shell effects in nuclear structure, quite well reproducible. 
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for example, by the model of independent particles in the potential 

well. 

2. If in the nuclear dynamics a chaotic component is dominant, it is neces¬ 
sary to expect, that the droplet model or Thomas-Fermi approximation 
will appear more useful. 

With such an approach, the elucidation of the mechanism of destruction 
of shell effects in the process of the R — C transition plays the key role. More 
appropriately the problem can be formulated in the following way: how does 
the shell dissolve with deviations from integrability or, conversely, how do 
incipient shell effects emerge as the system hrst begins to feel its proximity 
to an integrable situation? 

As has been mentioned above, the hnite motion of integrable Hamilto¬ 
nian system with iV-degrees of freedom, in general, is conditionally periodic, 
and the phase trajectories lie on iV-dimensional tori. In the action-angle 
variables (/, 9) the Hamiltonian is cyclic with respect to angle variables, 
H = Hq{I). Poincare called the main problem of dynamics the problem of 
the perturbation of conditionally periodic motion in the system dehned by 
the Hamiltonian 

H = HoiI)+eViI,e) 

where e is a small parameter. The essential step in the solution of this 
problem was the KAM theorem, asserting, that at turning on non-integrable 
perturbation, the majority of non-resonance tori, i.e. the tori for which 


N 

i=l 


dHo 

dli 


^0 


are conserved, distinguishing from the unperturbed cases only by small (to 
the extent of smallness of e) deformation. As was noted before, at dehnite 
conditions the KAM formalism allows us to remove the terms depending on 
the angle out of the Hamiltonian, using the convergent sequence of canonical 
transformations. When it succeeds, we hnd that perturbed motion lies on 
rather deformed tori, so that trajectories, generated by perturbed Hamilto¬ 
nian remain quasi-periodic. In other words, the KAM theorem reflects an 
important peculiarity of classical integrable systems: they conserve regular 
behavior even at rather strong non-integrable perturbation. In the problem 
of our interest, concerning the destruction of the shell structure of quantum 
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spectrum, the KAM theorem also will be able to play an important role. 
Considering the residual nucleon-nucleon interaction as non-integrable addi¬ 
tion to the self-consistent held, obtained, for example, in the Hartree-Fock 
approximation, we can try to connect the destruction of shell with the de¬ 
viation from integrability. The existence of shell structure at rather strong 
residual interaction (or at large deformation) can be due to the rigidity of 
KAM tori, contributing to survival of regular behavior. Such an assumption 
seems rather natural, especially if we take into account that the procedure of 
quasi-classical quantization [721 Eg itself as well as KAM theorem are based 
on convergence of the same sequence of canonical transformations. 

The aim of this section is to trace the evolution of the shell structure of 
QO Hamiltonian. In numerical calculations of this section it will be more 
convenient to use non-scaled version of the QO Hamiltonian (2.2) . The 
unperturbed Hamiltonian Hq we will understand as the Hamiltonian of two- 
dimensional harmonic oscillator with equal frequencies: 


Hq = H{a = 1, 6 = 0, c = 0). 

Its degenerate equidistant spectrum is well known. At switching on per¬ 
turbation, the degeneration disappears and the shell structure forms. The 
number of states, for example the states of F'-type (the numerical results, rep¬ 
resented below, are relative to the states of this type, while analogous results 
are observed for the states belonging to another irreducible representations 
of Ca^-group) for the given quantum number N is equal to 

No = ^{m + Mod{N,2)), 

where 

= h2N + Mod{N,3)), 

o 

and Mod{N, M) is the remainder of division of N by M. 

It is obvious, that eigenfunctions of exact Hamiltonian QO are no longer 
the eigenfunctions of operators N and L. Nevertheless, as numerical calcula¬ 
tions show, one can use the quantum numbers N and L for the classihcation 
of wave functions even at rather large nonlinearity. The measure of nonlinear¬ 
ity, at which such classification (i.e. the existence of shell structure) remains 
reasonable, is connected with quasi-crossing of neighboring levels. Quasi¬ 
crossing we shall understand as the approach of levels up to the distances of 
order of numerical calculations error. 
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Figure 7.7: a) Energy spectra for QO Hamiltonian (2.2) depending on values 
of parameter b for region C W = !?>. Points mark the quasi-crossings of 
levels. Bold solid line shows the dependence of classical critical energy E^r 
on parameter b. The dashed line shows the quasi-crossing region border. 
Arrow shows the point of quasi-crossing of levels with numbers fc = 40 and 
k = 41. b) Energy spectra for QO Hamiltonian (2.2) depending on values of 
parameter b for region C W = 3.9. There are no quasi-crossings of levels. 


The dependence of energy spectra of QO Hamiltonian on the parameter 
b for the values W = 3.9 and IF = 13 are represented in Fig j7.7[ As can be 
seen from Fig j7.7^ , for the PES with IF = 13 at the approach to the line of 
critical energy of the transition to chaos, defined according to the negative 
curvature criterion, the destruction of shell structure, which we understand 
as the beginning of multiple quasi-crossings, takes place. At the same time, 
for PES with W = 3.9 on Fig j7.7| 3p] quasi-crossings are absent even at larger 
nonlinearity than for W = 13. 

The destruction of shell structure can be traced for the wave functions, 
using the analogue of thermodynamic entropy Sk [751 [76], 

Sk = — ^2 lu 


^as we have shown in Section]^ the local instability is absent for that case and the 
classical motion is regular at all energies 
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The character of changes of entropy in regions -Ri and R 2 , corresponding to 
regular classical motion, correlates with the transition from shell to shell (see 
FigQ. Two effects are observed in the region C corresponding to chaotic 
classical motion. Firstly, quasi-periodical dependence of entropy on energy 
is violated, which testifies to the destruction of shell structure. Secondly, 
monotone growth of entropy is observed on average with going out on the 
plateau corresponding to the entropy of random sequence at the energies 
essentially exceeding the critical energy. 
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Figure 7.8: Entropy S as function of the state number. The straight lines 
connect points corresponding to shell classihcation on the basis of N. 



Chapter 8 

Wave Packet Dynamics 


The investigation of the time evolution of non-stationary states, i.e. of wave 
packets in quantum systems with classical analogues that admit chaotic be¬ 
havior, provides important information about QMCS. The localized quantum 
wave packet (QWP) is the closest analogue of the point in phase space, which 
describes the state of a classical system. However, such correspondence be¬ 
tween localized QWP and the classical point particle in the chaotic region is 
broken down for a very short time interval. Let us explain this fact by using 
Takahashi arguments m- We take two localized QWP 'ki(x) and T 2 (x) 
which are put in the chaotic region at an initial time to be slightly differ¬ 
ent from each other so that the difference between {Ti|x|Ti) and {\k 2 lxlT 2 ) 
(or (TiIpITi) and (^21^1^2)) is very small. We assume that in chaotic re¬ 
gion the localized QWP does not either extend in a certain time interval of 
the order l/\/h like that in regular region. Hence, following the Ehrenfest 
theorem, the packets d'i(x) and T2(x) move as classical particles and the 
distance between (\ki|x|Ti) and (\k 2 lxlT 2 ) (or (\ki|p|Ti) and (\k 2 lplT 2 )) is 
increasing exponentially in time. Let us consider the superposition 

T(x) = Ti(x) + T2(x), (8.1) 


which also becomes a localized QWP in the initial state. From that assump¬ 
tion, we can expect that the QWP (8.1) does not extend in a certain time 
interval of the order 1 /\/h . However, considering the exponential increment 
of the distance between (\ki|x|Ti) and (\k 2 lxlT 2 ), lk(x) (8.1) extends expo¬ 
nentially in the chaotic region and does not behave as a classical particle. 
This result is inconsistent with initial assumption and implies that in the 
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chaotic region localized QWP s (i.e. tl'i(x), \I' 2 (x) and d'(x)) extend expo¬ 
nentially as does the classical probability distribution in the first stage of 
time development. 

In order to describe such unusual behavior of QWP we should address 
to the concept of quantum-mechanical phase space. There are a few well- 
established schemes to introduce phase-space variables in quantum mechan¬ 
ics. In the present study we shall follow the procedure proposed by Weissman 
and Jortner Let us consider an initially localized QWP T characterized 
by coordinates q and momentum p, 


q= (T|q|T), p = (d'|p|d'). 


Now we introduce the coherent states |p,q), which in the coordinate x- 
representation are given by Gaussian wave packets 


N 


(x|p,q) = n (^^i) 

i=i 



( 8 . 2 ) 


In the study of a system of N coupled harmonic oscillators, it is convenient 
to choose for constants cxj the rms zero-point displacements 


= 


/ 


h 

rrijUj 


where rrij are masses and Uj are frequencies of the uncoupled oscillators. 
With this choice of aj, the coherent states a = |p, q) become the eigenstates 
of the harmonic oscillator annihilation operators 


aj\a) = aj\a) = 


(8.3) 


where 


= 


V2 




- + 'ft ft 


is a complex number. 

Using these coherent states, it is possible to introduce the following 
quantum-mechanical phase-space density 


P'i>(q,p) = 
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where ^ is any general wave packet. This quantum-mechanical coherent- 
state phase-space density may be regarded as a quantum analogue of classical 
phase-space density, since it satisfies an equation of motion where the leading 
term (when expanded in powers of h), corresponds to the classical Liuville 
equation. The stationary phase space densities are the following 


PE(q,p) = l(o|B>l' 


(8,4) 


where |q() is given in terms of Gaussian wave packet (8.2), while the eigenstate 
\E) is given by squares of the projections of the eigenstates on the coherent 
state, Eq.(8.3). Using the well-known expressions for scalar products {a\NL) 


we hnally obtain 


p^(q,p) = ^e-5(l“+l^+l“-lh^ 


C 


NL 


N,L 




( *ri_L I • *n_ *n-p\ 

[a^ a_ + ja_^_ «_ ), 


where 


with 


a± = -y={a2 T *tti) 

v2 


1 , , N + L N-L 

«i .2 = ri+ = —-—, n_ = —-—, j = ±1. 


The phase-space density ^^(q, p) is a function of the four real variables pi, qi 
and P 2 ,q 2 - We can get the contour maps of Pe{Pi,( 1 i',P 2 ,Q 2 ) in fhe {p 2 ,q 2 ) 
plane, taking = 0 and calculating pi from the relation 


H{pi,qi = 0,p2,q2) = E 

Quantum Poincare maps (QPM), obtained this way, constitute the quan¬ 
tum analogues of the classical Poincare maps and can be used for the search 
for QMCS both in the structure of wave functions of stationary states and 
in the dynamics of wave packets. 

Let us consider next the time evolution of a wave packet, which is initially 
in coherent state 

|4-(( = 0)) = |a) 

The time evolution of such initially coherent wave packet is given by 

= Y,\Ek){Ek\a)e-^^^^ 

k 
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The probability p{t) of finding the system in its initial state 
probability — is given by 

Pit) = \9a{t)\^ 

where Qait) is the overlap of ^(f) with the initial state 


the snrvival 
(8.5) 


Using the dehnition of stationary phase-space density (8.4), we can write 

3 


9ait) = 

k 


(8.6) 


Equation (8.6) implies that dynamics can be specihed by the spectrum of 
the initial coherent state |a). 

Weissman and Jortner ca have observed for Henon-Heiles Hamiltonian 
two limiting types of QWP dynamics of initially coherent Gaussian wave 
packets, which correspond to quasiperiodic time evolution and to rapid de¬ 
cay of the initial state population probability. Quasiperiodic time evolution 
is exhibited by wave packets initially located in a regular region, while rapid 
decay of the initial state population probability is revealed by those ones that 
are initially placed in irregular regions. An analogous situation takes place 
also in the multi-well potentials of our primary interest, in the energy region 
where the mixed state is realized. The quantum correlator p{t) for two-well 
potential reproduces well the characteristic properties of classical trajec¬ 
tories (see Fig j8.1[ ): quasi-periodic motion in the regular minimum manifests 
itself in oscillations of the correlator, and chaotic motion in another local 
minimum — in fast decay of correlations. 

Let us turn to the consideration of the dynamics of QWP in the PES with 
a few local minima (QO with W > 16) [2Z]. Transitions between different 
local minima can be divided into induced (the excitation energy exceeds the 
value of potential barrier) and tunnel transitions. The latter are subdivided 
into transitions from discrete spectrum into continuous spectrum (for exam¬ 
ple, a-decay, spontaneous division), and from discrete spectrum into discrete 
spectrum (for example, transitions between isomeric states). Up to now the 
process of tunneling across a multidimensional potential barrier, when ini¬ 
tial and hnal states are in discrete spectrum, has been the most complicated 
problem. 
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(b) 
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Figure 8.1: Typical trajectories for energies corresponding to the mixed state 
(a) and correlation functions (8.5) for equivalent Gaussian wave packets in 
chaotic (b) and regular (c) local minima of potential (2.4). 
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Most often the time evolution of wave packet is studied by two meth¬ 
ods: either by direct numerical integration of the Schrodinger time-dependent 
equation with corresponding initial condition 4/(r,t = 0), or by expansion 
of the packet 4/ (r, t) in the eigenfunctions of the stationary problem. The 
hrst one has some shortcomings, e.g., the complication relative to the inter¬ 
pretation of the obtained results and the necessity of separating the contri¬ 
butions from sub-barrier and tunnel transitions for the packet of an arbi¬ 
trary shape. These difficulties can be avoided if the sub-barrier part of the 
spectrum En < Uq {Uq is the height of the barrier) and the corresponding 
stationary wave functions t/’n(r) are known. Pure tunnel dynamics will take 
place for the packets representable in the form 

T(r,t) = En < Uq, 


Cn= ^{r,t = 0)^|Jn{r)dT. 


The probability p^{t) of hnding the particle at the moment of time t in 
certain local minimum R is 


p\t)= / |vI/(r,t)pdr = 5]C;ae 




^m(r)^n(r)dr 


R 

or in the two level approximation 
p^{t) =p^(0) -dCiCasin^ 
Let us introduce the value 


2h 


^jJl{r)^jJ2{r)dr. 


= max p^{t) 

Wt 

which is the maximum probability of hnding the particle in certain local 
minimum R, if initially it was localized to an arbitrary minimum. If the 
number of local minima is more than two, then of separateinterest is 

p^o — niinp'^“(t) 

Wt 

i.e., the minimum probability to hnd the wave packet in the minimum Rq 
corresponding to its initial localization. 
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Figure 8.2: Dependence of sub-barrier energy levels Ei for a double asym¬ 
metric one-dimensional rectangular well with infinite external walls (above) 
and on d (below). Width of the well equals 3, width of the barrier — 1, 
barrier height — 2. 


Intuitively, we may suggest that 1 if the initial minimum is lo¬ 

cal, and the hnal one is absolute. However, the results of pg obtained for 
the simplest one-dimensional models (asymmetric double wells of different 
shapes) are inconsistent with the intuitive expectations. The probability of 
tunneling from the local minimum to the absolute one depends resonantly 
on the potential parameters. Fig, 8^ gives the dependence of on the well 
depth displacement d. It can be seen that at an arbitrary asymmetry -C 1. 

The resonant behavior of p^ becomes more clear if one considers the 
spatial structure of the sub-barrier wave functions. For a sufficiently wide 
barrier in the case of an arbitrary asymmetry, the sub-barrier wave functions 
are largely localized in separate minima. The delocalization takes place only 
in the vicinity of the level quasi-crossing. The degree of this delocalization 
directly depends on the distance between the interacting levels. Obviously, 
the QWP, in which the components localized in the certain minimum are 
dominating, cannot tunnel effectively to the neighboring minimum. It is 
precisely this reason that explains the stringent correlation between the p^ 
minima and the level quasi-crossing (see Fig,8.2). 
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Now the question arises if a similar correlation between the level quasi¬ 
crossings, the delocalization of wave functions and resonant tunneling persists 
in the two-dimensional case. To give an answer to this question let us turn 
to QO Hamiltonian (2.2) in the region W > 16. Recall that three {C^y- 
symmetric) identical additional minima appear at hh > 16 apart from the 
central minimum. The central minimum exceeds the lateral ones in depth 
in the region 16 < IT < 18. At IT > 18, the central minimum becomes 
the local one. In this region of parameters the procedure of diagonalization 
of QO Hamiltonian in oscillator basis becomes essentially complicated. It is 
connected with the fact that the basis of Hamiltonian, the potential of which 
has the unique minimum, is used for the diagonalization of Hamiltonian 
with complex topology of the PES. In addition to the large dimension of the 
basis, it is necessary for the basis wave functions to have sufficient value in 
the region of lateral minima. It has been achieved by the optimization of 
frequency ojq of the oscillator basis. For the values of parameters used in 
the calculations (IT = 17.8, h = 0.17) the value of the optimal oscillator 
frequency was lcq = 0.2. Fig j8.3| gives lower eigenvalues of E and Ai types 
depending on the basis dimension. We can see that it is possible to get 
saturation in the basis at the dimension of sub-matrices N ~ 10^ even for 
low states localized in the lateral minima (dashed lines). 

We can also use the above introduced analog of thermodynamic entropy 
Sk to estimate the degree of saturation in the basis. Fig,^ gives Sk values 
for the states of Ai type for the dimension of sub-matrices 408 and 690. 
Increasing of the basis dimension does not lead to signihcant changes of Sk 
values for the states with energy up to the saddle one Es- 

The states localized in the central or in the lateral minima have essentially 
different distributivity of the coefficients C^, {i = NLj} (see Fig j8.5 r,b) and 
thus different entropies: states localized in the central minimum have less 
entropy. In the neighborhood of the points of level quasi-crossings, delocal¬ 
ization of wave functions corresponding to these levels takes place; these wave 
functions have close distributivity of the coefficients Cf (see Fig j8.5| 3,d). 

Fig,8T shows the sub-barrier part of the energy spectrum obtained by 
the diagonalization. As is easy to see, the tunneling of the wave packet, 
composed of the sub-barrier wave functions, can be described in a two-level 
approximation. Indeed, there are approximately 10 level quasi-crossings of 
Ai and F-types, where the nonlinearity parameter b changes are of the order 
of 10“^, while the effective half-width of the overlap integral in (3.66) is about 
10“® (see Fig|8.7[). Hence, all non-diagonal elements will be close to zero with 
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b 


Figure 8.3: Dependence of the energy spectrum of E-tjpe (a) at 6 = 

0.173798 and Ai-type (b) at 6 = 0.17 on the dimension of sub-matrices N 
for W = 17.8. Dashed lines correspond to the states localized in peripheral 
minima. 
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Figure 8.4: Dependence of Sk on the state number k (Ai-type) for the 
dimension of sub-matrices 408 (dark points) and 690 (crosses) at W = 17.8 
and 6 = 0.17. Es corresponds to the saddle energy, and Ecr — to the 
classical critical energy of the regularity-chaos transition. 




Regular and Chaotic Classical and Quantum Dynamics 


147 
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Figure 8.5: Distribution of coefficients C^ by number of basis state i = 
{NLj) {E - type) for localized states with /c = 3 (a) and /c = 4 (b) at 
b = 0.17379 {W = 17.8) and for the delocalized states with k = 3 (c) and 
A; = 4 (d) in the point of quasi-crossing b = 0.1737924 {W = 17.8). 
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Figure 8.6: Dependence of energy levels Ei of the QO Hamiltonian (2.2) on 
the parameter b for PF = 17.8: a) — spectrum of i7-type, b) — spectrum of 
Hi-type. 


the overwhelming probability in the matrix 

ttmn = j V’m(r)^/’n(r)dr 

R 

of an arbitrary nonlinearity parameter (e.g., b). Two essential non-zero out- 
of-diagonal matrix elements (two-level approximation) appear only in the 
vicinity of quasi-crossings. The probability of double quasi-crossings at a 
hxed nonlinearity parameter is almost excluded. This probability is by two 
or three orders lower than that of rather rare ~ 10“^ single quasi-crossings. 

So now we can give the answer to the question we have posed above. The 
stringiest correlations between level quasi-crossings, delocalization of wave 
functions and resonant tunneling across the potential barrier take place in 
the two-dimensional case (and, most likely, in the multi-dimensional one too). 

The existence of the mixed state for many-well potentials must be clearly 
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Figure 8.7: a) - Quasi-crossing of energy levels of i7-type with A: = 3,4 for QO 
Hamiltonian (2.2); b - localization of wave functions of states with fc = 3,4 
in the central well at different 6; c - dependence p^{h). 
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manifest in the dynamics of quantum wave packets. The pre-exponential 
factor of the tunnel amplitude depends on the type of classical motion and 
consequently we expect to observe asymmetry of the effective barrier pene¬ 
tration in the mixed state. This purely classical effect may be observed only 
if the uncertainty in level energy is comparable with the average distance 
between levels, and the system does not yet ’’feels” that the spectrum is dis¬ 
crete. This determines the time scale for which observation of the effect is 
possible. It is the same time scale on which transition from classical diffusive 
evolution with linearly increasing energy to the quasi-periodical quantum one 
is observed [5U] . 


8.1 Dynamical tunneling 

Transition from integrability to chaos in classical dynamics not only signih- 
cantly modifies the tunneling process, but moreover it leads to the appearance 
of principally new scenarios of the tunneling effect. The simplest example 
of such a scenario is the so-called dynamical tunneling [81]. Dynamical tun¬ 
neling appears in systems where the phase space contains such regions that 
transition between them is forbidden in classical mechanics, but this restric¬ 
tion is not due to the potential barrier. Of course such an effect is possible in 
systems with two or more degrees of freedom, where additional to the energy 
integrals of motion are sources of the corresponding restrictions. This new 
type of tunneling is more complicated than the traditional (barrier) tunnel¬ 
ing. The complication is connected with the fact that simple consideration of 
the potential energy surface does not reveal the conditions of the restriction. 
Instead of static potential surfaces one should consider dynamical behavior 
of the trajectories. 

In order to understand the origin of dynamical tunneling in a bounded 
system, let us recall the semiclassical solution of the one-dimensional sym¬ 
metric double well problem. If we quantize the system, considering each well 
separately, we obtain a spectrum composed from exactly degenerate doublets. 
Only by taking into account the interaction between the wells via overlap of 
exponentially small sub-barrier tails of the wave functions, we can obtain 
the correct result — almost degenerate pairs of levels with the well-known 
splitting. The magnitude of the latter determines the tunneling rate. 

A situation similar to those in the one-dimensional symmetric double well 
can exist also in a multi-dimensional potential in the absence of energy bar- 
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rier. Let us consider [821 [83| l8^ a dynamical system with reflective symmetry 
of the phase space T. Let there be two non-connected regions Ai and A 2 in 
the phase space, each of them being invariant with respect to classical dynam¬ 
ics, mapping one onto another with the symmetry transformation A 2 = TAi. 
Let the classical motion in 2 be regular, i.e. those regions represent the 
stability islands embedded in the chaotic sea. In that case as we have seen 
above the regions ^41^2 represent invariant tori in classical phase space. We 
make an additional assumption that in the semiclassical limit there is a set 
of states '0(q) mostly localized in the regions Ai and A 2 . Using the standard 
procedure we can quantize motion in each region independently and con¬ 
struct degenerate wave functions cind sometimes 

called quasi-modes. Taking into account the interaction between the regions, 
the quasi-modes one-dimensional case, must be replaced 

by their combinations — symmetric and antisymmetric: 

Energy degeneration between the functions is removed due to the tunnel¬ 
ing processes. However, in contrast to the one-dimensional case, the invariant 
tori Hi 2 are not always separated by the energy barrier in the configuration 
space. Transitions Hi o H 2 may be classically forbidden by conservation of 
the constant of motion distinct of the energy. 

Like in the one-dimensional case, a wave packed constructed from linear 
combinations of the wave functions 


='^^■'■^(q) ±^(q) 


(8.7) 


will mimic a classical particle, initially localized in one of the regions Hi_ 2 . 


The wave packet (8.7) will make the transitions between regions with the 
frequency determined by the energy splitting of states and '0*'“^(q). 

A natural question arises: if in the phase of the system there are two sym¬ 
metric stability islands separated by chaotic sea, then how that sea affects 
the dynamical tunneling between the islands? Semiclassical reasoning points 
out that the probability distribution associated with the wave function of 
the quantized torus decays exponentially outside the torus. The small over¬ 
lap in the classically forbidden region of the decaying distributions, centered 
on the two congruent quantized tori, will lead to the tunneling splitting. If 
there is no chaotic region between the tori, then the overlap will be very 
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small. However if a chaotic region exists between the tori, then the wave 
fnnctions corresponding to the tori overlap initially with the chaotic states. 
Due to the ergodic nature of chaotic wave functions (uniform distribution of 
the probability density), the connection between the two tori appears to be 
more efficient that in the case of a regular intermediate state. Therefore we 
can expect that the tunneling will be assisted by the presence of the chaotic 
region. 

All these qualitative considerations are conhrmed by the analysis [82] of 
the splitting magnitude in the tunneling doublet as a function of the pa¬ 
rameter a — the measure of chaoticity of the system. The splitting under¬ 
goes considerable fluctuations, which are connected with quasi-crossing of 
the ’’outside” chaotic level and the considered regular tunneling doublet. It 
is important to note that because there is no dynamical fragmentation of the 
chaotic region on separated symmetric parts as happens in the regular part 
of the phase space, we have no reason to expect the appearance of chaotic 
doublets. So every chaotic state has hxed parity. 

Peculiarities of such a spectral region can be described in the three-level 
model [83] dehned by the Hamiltonian 

/E+e 0 0 \ 

H = \ 0 E — e V I (8.8) 

\ 0 n E^ } 


Here EEe are energies of the regular quasi-doublet, composed of the symmet¬ 
ric and antisymmetric combinations of the corresponding quasi-modes. The 
model assumes that the even chaotic state \C) with energy E^ is connected 
with the |-|-) state with energy E — ehy the coupling constant v. In practice 
one of the constants u oi e dominates. If £ dominates then we turn back 
to the above considered two-Ievel model. Therefore let us assume that the 
dominating constant is the coupling between the \C) and |-|-) states and let 
us take £ = 0. In that case the doublet splitting AE will be governed by the 
energy shift £'+ due to the coupling v. Diagonalization of the Hamiltonian 
(8.8) leads to 


AE = 


e-ec 

\v\ 


E-E^ 

E-E^ 


It means that with the variation of a in the tunneling doublets splitting we 
will observe peaks with magnitude |z/| when E^{a) crosses E^{a). 

Now we draw some conclusions. As distinct from the integrable case, 
systems with dynamics of mixed type, where the phase space contains both 
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regular and chaotic regions, demonstrate a new tunneling mechanism called 
chaos assisted tunneling [82]. The doublets splitting that determines the 
tunneling rate in two-level approximation in systems of mixed type, is as a 
rule several orders of magnitude higher than for similar integrable systems. 
In contrast to direct processes when a particle tunnels directly from one state 
to another, chaos assisted tunneling corresponds to the following three-step 
process: 


1. tunneling from a periodic orbit to the closest point in the chaotic sea; 

2. classical propagation in the chaotic region of the phase space towards 
the vicinity of another periodic orbit; 

3. tunneling from the chaotic sea onto the latter periodic orbit. 


In other words doublet splitting due to reflective symmetry occurs not di¬ 
rectly but through the compound-process of wave function destruction, piece 
by piece, near one regular region, and then chaotic transport through the 
chaotic sea towards the neighboring symmetry conjugated regular region, 
and at last restoration of the initial state image. Schematically this pro¬ 
cess is presented on Fig j8.8 Let us note that chaos assisted tunneling is a 
formally higher perturbative order process than direct tunneling. But the 
corresponding matrix elements for chaos assisted tunneling are much larger 
than for the direct process. Intuitively it can be pictured as the following: 
as the main part of the distance (through the chaotic sea) represents a clas¬ 
sically allowed transition, then we can expect that those indirect trajectories 
will make a larger contribution to the tunneling flow than the direct ones. In 
the case of direct tunneling all sub-barrier trajectory represents a classically 
forbidden process. 

Let us now consider the dynamics of Gaussian wave packets in QO po¬ 
tential (2.3) with W = 18, when PES have four local minima — the central 
one with chaotic motion and three symmetric peripheral with regular motion 
for sub-barrier energies (see hg,8.9r). 

At considerably super-barrier energies of motion, chaotic dynamics is ob¬ 
served in all local minima; however in the peripheral one, a large stability 
island survives even for E = 2Es (hg,8.9D), which correspond to quasiperi- 


odic trajectories localized in those minima (hg,8.9 


With the spectral method it is easy to perform a numerical simulation 
of the dynamical tunneling process for super-barrier motion in QO potential 
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Figure 8.8: A schematic representation of direct and chaos assisted tunneling 




Figure 8.9: Regular and chaotic trajectories for E < Es (a), PSS (b) and 
regular quasiperiodic trajectory (c) for E = 2Es in the QO potential (2.3). 
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(2.3). In order to do that we prepare the initial state ipo{x,y) in the form 
of a Gaussian minimum uncertainty wave packet (8.2) with such parameters 
{xo,yo,Pxo,Pyo) that correspond to classical initial conditions for a trapped 
regular trajectory (£g,8.9::). For some time the quantum dynamics of such a 


wave packet will imitate classical motion on the semiclassical trajectory. Af¬ 
ter a sufficiently long period of time, the dynamical tunneling effect starts to 
be visible — the wave packet leaves the local minimum of its initial localiza¬ 
tion and appears in two others — a process forbidden in classical mechanics. 
After the same time period the wave packet more or less restores in its initial 
states (£g,8.10). 

The oscillating nature of dynamical tunneling is clearly described in terms 
of the autocorrelation function (5.23) (£g,8.11 on the left). Its Fourier trans¬ 
form (hg,8.11 on the right) demonstrates three dominant states in the wave 
packet composition, that form the tunneling triplet (the right maximum on 
the P{E) diagram is actually a very tiny tunneling doublet). Figure 8.12 
shows the probability distribution \'ijj{x, y)\‘^ for all three members of the tun¬ 
neling triplet in order of energy increasing: the intermediate state is purely 
regular and forms a tunneling doublet with the state of highest energy which 
is almost regular with a small admixture of chaotic modes. The lowest lying 
state is predominately chaotic and overlaps with only one of the tunneling 
doublet partners (the higher one). Repulsion of those states considerably 
broadens the level splitting in the tunneling doublet and therefore increases 
the tunneling rate substantially, which is the essence of the chaos assisted 
tunneling effect. 
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Figure 8.10: Temporal evolution of gaussian wave packet (8.2) under dy¬ 
namical tunneling in the QO potential (2.3). 



Figure 8.11: The auto correlation function (5.23) (left) and its Fourier 
transform (right) for Gaussian wave packet (8.2) tunneling dynamics in the 
QO potential (2.3). 
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Figure 8.12: Level lines for |'^(a:,|/)p distribution for the tunneling triplet 
wave functions (the energy grows from left to right). 












































Chapter 9 

The mixed state and concrete 
physical effects 

9.1 Decay of mixed states 

The escape of trajectories (particles) from localized regions of phase or config¬ 
uration space has been an important topic in dynamics, because it describes 
the decay phenomena of metastable states in many branches of physics: 
chemical and nuclear reactions, atomic ionization, nuclear fusion and so on. 
This problem has a rich history. Almost a century ago, Sabine |86] considered 
the decay of sound in concert halls. Legrand and Sornette |SZ] have shown 
that this problem is equivalent to the escape one: a small opening of width A 
for escape must be identified with f a(S)ds, where a(s) is the absorption co¬ 
efficient at position S of the container (billiard) boundary, q;(s) = 1 over the 
width of window and q;(s) = 0 elsewhere. Szepfalusy and Tel [88] connected 
the escape problem with the problem of chaotic scattering. 

Exponential decay is a common property expected in strongly chaotic 
classical systems [83 EHl E] • Let us consider as an example (82] point par¬ 
ticles in a rectangular box bouncing elastically off the walls. We allow our 
system to decay by providing a small window in one of the box walls through 
which particles are allowed to escape. As is well known, the motion of par¬ 
ticle in a rectangular billiard is regular: two independent integrals of motion 
are the absolute values of momentum projection on the billiard walls. The 
trajectories of particles become chaotic if a circular scattering center is placed 
somewhere inside the box. 
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For the chaotic case, simple consideration leads to exponential decay. The 
number of particles leaving per time interval is given by 

^ = Ap{t) j (fp pe„ = -2Ap{t)p^5p (9.1) 

Here p is the absolute value of momentum, e„ is a unit vector normal to the 
opening in the surface, and integration in momentum space is taken over a 
circular ring with radius p and inhnitesimal width 5p. Function p{t) is the 
phase space density, which in the ergodic motion is only a function of time. 
In our case 

Nit) 


pit) = 


(9.2) 


2TTp5pAc 

where A^ is total coordinate space area available. Inserting (9.2) into ( 9.1| ) 
yields 

vA 

N{t) = iV(0)e-“*; a = ^ (9.3) 


TlAr 


Analytically calculated decay constant a is in a good agreement with the 
graphically extracted value. 

Exponential law at extremely long times turns into the power law typi¬ 
cal for decay of regular systems. One possible mechanism for generation of 
power tails is the effect of ’’sticking” of the chaotic orbits to outer boundaries 
of stability islands [92] , or a very similar effect connected with the existence 
of marginally stable periodic or ’’bouncing ball orbits”. Although some qual¬ 
itative models, which show how the algebraic tail emerges, were introduced 
in [HS], no critical conditions for the distinct decay laws were formulated in 
terms of the billiard geometrical constrains. Experimental escape of cold 
atoms from a laser trap of billiard type with a hole was studied in [9H 1^ . 

Transition from the billiard escape problem to escape from potential wells 
from one hand substantially broadens the number of possible applications, 
and on the other signihcantly complicates the problem. Of course, the one- 
well case is a simpler one. Zhao and Du (HS] reported a study on the escape 
rates near the threshold of Henon-Heiles potential (3.3). Simulations per¬ 


formed by the authors show that the escape from a Henon-Heiles system at 
energy slightly exceeding the saddle one, follows an exponential law similar 
to chaotic billiard systems. They derived an analytic formula for the escape 
rate as a function of energy. The derivation is based on the fact that the 
phase space considered potential (as well as the billiards) is practically ho¬ 
mogeneous near the saddle points. It should be noted that in such a case all 
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trajectories with energy higher than the saddle one leave the potential well 
in hnite time. The only problem to solve is to determine the probability of 
a particle to escape from the well in unit time interval. 

We intend to study the particles escape from the local minima in the 
case when the phase space contains macroscopically signihcant components 
of regular as well as of chaotic type EZI- As we have seen, such a possibility 
is realized in potentials with several local minima in the form of a mixed 
state. The problem of particles escape from a potential well in the presence 
of the mixed state has an essential distinction from the one considered in 
|nS]- Detailed analysis of the Poincare sections for the potentials and 
QO, made in chapter shows that in the ’’regular” well at energy higher 
than the saddle, chaotic trajectories appear and their measure increases with 
energy growth. At the same time the regular trajectories remain mostly 
trapped in the well. The phase space of such a system contains macroscopic 
regions of both the regular and chaotic motion. Therefore the ensemble of 
particles initially situated in the well, divides on two components. So the 
ultimate problem is to determine the probability to escape per unit time for 
the chaotic component, and the relative measure of non-escaping trajectories 
trapped in the well. 

Let us address again the example to above considered in details potentials 
of quadrupole oscillations (2.3) and umbilic catastrophe (2.4). In the hrst 
as well as in the second case, in some local minima up to the saddle energies 
the motion remains absolutely regular. Moreover, at energies signihcantly 
higher than the saddle energy (see £g 9.1) the phase space structure preserves 
division on chaotic and regular components. The latter is localized in the 
part of conhguration space which corresponds to regular motion at energies 
below than saddle energies. 

This phase space structure leads to the fact that escape from such local 
minima has all the above mentioned properties of decay of chaotic systems, 
and also a diversity of principally new features, representing an interesting 
topic for conceptual understanding of chaotic dynamics, and for application 
as well. We are interested only in the ”£rst passage” effects, leaving aside the 
problem of dynamical equilibrium setup for the finite motion (for example, 
in QO potential). It is important to stress that though we study the process 
of escape from a concrete local minimum, the over-barrier in the case of a 
mixed state has a specihc memory: general phase space structure at super¬ 
saddle energies is determined by the characteristics of motion in all other 
local minima. 
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Figure 9.1: Poincare surfaces of section ?/ = 0 in the {x,px) plane for the 
potentials (a) and QO (b) for E = 2Es 

We carried out numerical simulation and analytical estimates EH of tra¬ 
jectories escape in potentials and QO through the hole over the saddle 
point. Results of the escape problem for systems with multi-component phase 
volume (regular and chaotic components) essentially depends on choice of ini¬ 
tial ensembles for dynamical variables. Fig j9.2| presents the normalized par¬ 
ticle number N{t)/N{t = 0) for 10® initial conditions, uniformly distributed 
inside the right minimum in the potential (a: > 0), and peripheral mini¬ 
mum in the QO potential {x > 1/12), together with the typical trajectories 
and Poincare sections. The results for different potentials are evidently sim¬ 
ilar and have such characteristic features: 

• At times t —>■ oo decay law saturates at 

N{t^oo)= 

where is equal to the relative phase volume of ’’never-escaping” 
trajectories, which represents the regular trajectories completely local¬ 
ized inside the considered minimum. All such trajectories, therefore, 
have inhnite escaping times. 

• For t > t{E) the decay law has exponential form 

N(t)/No = p<"'> + 


(9.4) 
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Figure 9.2; Decay law for mixed states in the (a) and QO (b) poten¬ 
tials. Solid lines — numerical simulation for E/Eg = 1.1,1.5, 2.0; dotted 
and dashed lines — theoretically obtained exponential and linear decay laws 
respectively, zoomed on the inset hgure in the lower left corners. Other inset 
hgures represent the typical trajectories and Poincare sections for the three 
different types of initial conditions: linearly escaping, exponentially escaping 
and non-escaping. Circles show the joining points between the linear and the 
exponential decay laws at critical time t = r. 
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where represents the relative number of exponentially escaping par¬ 
ticles. 


t < t{E) For the decay law is linear 

N{t)/No = 1-a^‘h. 


(9.5) 


We should stress, that (9.5) is in no way a linear approximation of (9.4) 


for small f: in general e ^ 1 anda^^^p^^^e ^ Instead, 


from the condition of smooth joining of curves (9.5) and (9.4) in the transition 
point t = T, we obtain 

p(e) = 1 — pG'^) — pih 

where is the relative number of linearly escaping particles. More¬ 


over, already on time scales t < t{E) linear decay law (9.5) is apparently 
different from its exponential analogue 

N{t)/N{t = 0) = p("") + (1 - 


(see the inset on £g,9.2). 


As we can see from the inset Poincare section on £g.9.2 both the chaotic 


and regular trajectories contribute to linear escaping regime (9.5), because 


for sufficiently small times t < r chaotic and regular motions are not yet 
distinguishable. Up to transient time f = r all quasi-one-dimensional regular 
trajectories, oriented along the x-axis, have already escaped and for f > r 


the escape of the remaining chaotic particles follows exponential law (9.4). 


The particles escaping the last show already mentioned sticking phenomena 


(see £g,9.2). 


The transient time t{E) in fact coincides with the passage time of the 
longest one-dimensional path from the opening to the opposite wall of the 
potential well and back (see fig. 9.2[ ). For the potentials and QO corre¬ 
sponding theoretical estimates read (we assumed m = 1) 


Vai+Ue) 

Tn.{E)^2 J ifj = 4A- 

0 

Tqo(E) ^ 12 (fy K 





Ve 


= 6V2td,(|:) 




















Regular and Chaotic Classical and Quantum Dynamics 


163 


where K{k) is the complete elliptic integral of the first kind and Es = 1/12^ 
is the saddle energy in the QO potential for W = 18 (for the potential 
Es = 1). 

Theoretical estimates for the escape rate were obtained by averaging the 
escape probability over the opening PB] : 


a{E) = p{E) J dy J d9\p\cos9, 

x=xs -f 

where X 5 is the coordinate of saddle point and p{E) is the normalized particle 
density: 

^ 271 A{E)' 

where A{E) denotes area of the classically allowed region inside the well: 


A{E) = 


dxdyQ {E — U{x, y)). 


x>xs 


Such p{E) corresponds to uniform distribution of particles on the energy 
surface H{p,q) = E. 

For the potentials and QO the explicit formulae are the following: 


a/2(1+\/F) 


Ary{E)=2 


dx\ 


X + 2a 


i+J— 

Eg 


^qo{E) - — 




0 


where ^ = \fEs{x — xs)- 

Finally, the general expression for the escape rate is 


a{E) 


(P) 

7rA{E)' 
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For our case (p) = f dy\p\. In the case of billiards with small opening 


x=xs 


p = const ^ (p) = pA, and we recover expression (9.3). For the potentials 
and QO we get the results in closed form: 




E-1 


2y/aA]:,^{E) 


^ (Z^) - (Z^)} 

where e = E — Es + 1/256, K{k) and E{k) are the complete elliptic integral 
of the hrst and second kind respectively. 

In order to obtain ad\E) we correct A{E) subtracting the relative phase 
space occupied by the non-escaping particles 


= 4(B)(1 - ^ aV\E) = 


Fig |9.2| demonstrates good agreement between our theoretical and numerical 
results for a wide energy range. 

The fraction of the non-escaping particles coincides with the relative 
phase space volume of trajectories, localized in the regular minimum, which 
may be well estimated by the relative area of the stability island on the 


Poincare section (£g,9.2). Calculation of the relative area of regular island in 
the Poincare section was performed by the following scheme. By numerical 
integration of equations of motion the island boundary was determined and 
then the interior area was calculated. Further the obtained area was divided 
on the entire area of classically allowed motion, dehned by the conditions x > 
0 and p^ > 0. While the phase volume itself is 4-dimensional, the stability 
island in Poincare section is 2-dimensional, and so we cannot expect absolute 
coincidence of the corresponding measures. However, the calculations show 


very close correspondence between them (see fig,9.3). Therefore, numerical 
analysis of Poincare sections together with our theoretical results gives all 
the information necessary to predict the escape dynamics by an independent 
method. 

In summary, we have considered classical escape from separated local 
minima in two representative 2D multi-well potentials, realizing the mixed 
state. We have found that escape from regular minima contains a number of 
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Figure 9.3: Correlations between relative area of the stability island 

and non-escaping particles ratio empty squares - filled squares 

- p^^^\e), empty circles - Pqq{E), filled circles - Pqq{E). 
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new features. The most important among them are the following: i) decay 
law saturates at long time ranges; ii) on small time scales there exists a linear 
segment, which is not connected with linear approximation to the exponential 
decay law observed in chaotic systems with homogeneous phase space. The 
fraction of particles, remaining in the well is determined by the relative phase 
volume of the regular component, which in its turn monotonically decreases 
with growth of energy. It was shown that the linear segment of the decay 
law is generated by the quasi-one-dimensional trajectories, oriented perpen¬ 
dicular to the opening, and the transient time of the linear-to-exponential 
regime lies in perfect agreement with the analytical estimates. 

We should note that we devote our main attention to escape from the 
regular local minima, because the specihcs of the mixed state manifests only 
in them. However let us recall that in the case of mixed state the phase space 
structure at super-saddle energies is determined by dynamical characteristics 
in different local minima of the whole potential energy surface. 

The above mentioned peculiarities of the escape problem may have found 
practical application for extraction of the required particle number from 
atomic traps. Changing the energy of particles trapped inside the regu¬ 
lar minimum, we can extract from the trap any required number of parti¬ 
cles. The problem of particle energy changing in the potential well may be 
solved by the introduction of small dissipation. The obtained results may 
be of interest also for the description of induced nuclear hssion in the case 
of double-humped hssion barrier. Revealed peculiarities must manifest also 
in over-barrier dynamics of wave packets, initially localized in the regular 
minima. 


9.2 Chaos assisted tunneling in superdeformed 
nuclei 

Since its discovery the atomic nucleus has been constantly used for verifying 
new physical ideas such as superhuids, superconductivity, supersymmetry 
and dynamical chaos. The aim of this section is to demonstrate the possibility 
of observing the above mentioned features of chaotic multi-well dynamics in 
concrete physical effect — the decay of the superdeformed states of atomic 
nuclei. In contrast to the classical decay of the mixed state considered above, 
this effect is purely quantum, but nevertheless dependent on the structure of 
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the classical phase space. 

The superdeformed (SD) nucleus has an ellipsoid shape with axes ra¬ 
tio a/h ~ 2. Such extreme nuclei shapes are the result of the quantum 
nature of particles that compose the nucleus and £11 the discrete energy- 
levels (shells). In large deformations, the energy interval (gap) between 
filled and next free shells could be very large. In combination with other 
effects this energy gap could provide stability of the nucleus in such large 
deformations. Superdeformed nuclei were found more than twenty years 
ago [HB], and today they have been observed in several mass regions around 
A = 20, 40,80,130,150,165,190 and 240 


Study of SD nuclei offers a new way of understanding nuclear structure. 
This is a good illustration of the above mentioned fact that nuclei are a 
laboratory for research of the general physical effects such as tunneling, chaos 
and phase transitions. On the one hand, nuclei provide interesting data for 
the study of these general phenomena. On the other hand, we obtain new 
information about nuclear structure. 

SD nuclei are produced in reactions with heavy ions. Initially accelerated 
heavy ions, when colliding with nuclei of target material, produce highly ex¬ 
cited and fast rotating compound nuclei. These nuclei release part of the 
excitation energy by the emission of light particles (neutrons, protons, alpha 
particles) and photons. In all observations of fast rotating SD nuclei, rota¬ 
tion breaks off (already at low momentum), when nucleus suddenly changes 
shape and decays to the state that corresponds to lower deformations. Three 
stages of the transition from SD state to the normal deformed (ND) one 
are presented on the fig raHODl: feeding of superdeformed bands, ordered 
rotation and decay from superdeformed to normal states. 

The lowest states in SD local minimum correspond to high excitations 
in the main minimum. As a result, the lowest SD states (cold states) are 
characterized by good quantum numbers and symmetries, while their decays 
to the main minimum are controlled by strict rules of selection. In contrast, 
ND states at the same energy could correspond to chaotic motion in the 
semi-classical limit. Statistical analysis confirms this idea: energy spectrum 
in the region Eext ~ 0 — 7Mev demonstrates all the signs of chaos, while 
yrast region is regular nnn. This means that in the considered region of 
nuclei excitation for which deformation potential allows for the existence 
of second minimum, there is realized a mixed state, considered in detail 
through the present report. In that case a natural question arises: can the 
chaos assisted tunneling mechanism be used for description of the transition 
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Figure 9.4: The three stages in the ’’life” of a superdeformed nucleus together 
with some characteristic properties unni. 
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from the SD state to the ND one? We should stress that we can consider 
dynamical tunneling because the transition takes place in multi-dimensional 
space: nuclear shape is characterized by at least two parameters. 

Aberg |1U2] estimated chaos assisted tunneling in perturbation theory. In 
the limit of no mixing between the ND states (completely regular ND system) 
the wave function of the SD state mixes with the doorway state only. In this 
case the tunneling probability is given by 


rjnregular 



2 


where V) is tunneling coupling between SD and doorway states, and 2Ai? is 
the energy distance between the doorway states. 

In the other limit situation (quantum chaos) the tunneling strength is 
spread out over all ND states. Because tunneling probability in the chaotic 
case is 

rpchaotic 


Consequently 

^chaotic 

-~ TV ~ --—- 

'T'regular ^ , 

Hdoorway 

For pdoorway ~ 1 Mev~^, and relative excitation energy between the SD 
and ND state E = 3 — 5 Mev Pnd/ Pdoorway ~ 10^ — 10®. Therefore we expect 
the tunneling probability to be enhanced by 10^ — 10® times if the ND states 
are chaotic. In other words, the tunneling process connected to the decay 
out of superdeformed states is strongly enhanced by chaotic properties of the 
ND states. 
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Chapter 10 
Conclusions 


Autonomous Hamiltonian systems with many local minima describe a variety 
of physical processes, including chemical nuclear reactions, nuclear hssion, 
phase transitions etc. Nevertheless, the potentials of non-trivial topology 
still remain far aside of attention, from both the technical and ideological 
difficulties. The present work is one of the hrst studies of regular and chaotic 
classical and quantum dynamics in multi-well potentials. 

The most general type of classical motion in potential with two and more 
local minima represents the so-called mixed state, which was hrst observed for 
simple polynomial 2D potentials. We proposed quantum mechanical treat¬ 
ment of the mixed state, and demonstrate signatures of quantum chaos in 
the energy levels statistical properties, stationary wave function structure and 
Gaussian wave packets dynamics. In particular, we have presented a new ap¬ 
proach to the investigation of QMCS in wave function structure, which can 
be realized in potentials with two or more local minima. The main advantage 
of the proposed approach is the possibility to detect QMCS in comparison 
not different wave functions, but different parts of the same wave functions. 
Efficiency of the approach was demonstrated for the deformational potential 
of surface quadrupole oscillations of nuclei and lower umbilic catastrophe D^. 

In this work we study the energy spectra statistical properties in the 
case of the mixed state. Spectral series in the energy ranges, where the 
mixed state is realized, open new possibilities for studies of intermediate 
statistics. At these energies chaotic and regular components are separated 
not in phase space (as usually the case), but in the conhguration space. A 
priori the nearest neighbors spacing distribution function in the mixed state 
does not obviously reduce to the properly weighted superposition of Poisson 
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and Wigner distributions. In the mixed state we deal not with statistics 
of a mixture of two level systems with different nearest neighbors spacing 
distribution functions, but with statistics of a level system, where each level 
does not belong to dehnite statistics. Statistical properties of such systems 
still wait investigation, however such systems correspond to the common case 
situation. 

The solution of Schrodinger equation for two-dimensional potential sys¬ 
tems of general form imposes high requirements on efficiency of numerical 
methods used. The matrix diagonalization method is efficient for one-well 
potentials. However, this numerical procedure becomes less attractive at the 
transition to multi-well potentials. We have shown that, in this case an at¬ 
tractive alternative to the matrix diagonalization may become the spectral 
method. The advantages of the spectral method compared to the diagonal¬ 
ization method are demonstrated in the analysis of quantum chaos problems 
for two-dimensional potentials with complicated geometry. 

The report presents also some concrete realizations of physical effects, 
connected with the mixed state. Having considered classical escape from 
2D multi-well Hamiltonian system, realizing the mixed state, we show that 
escape from local minima has a diversity of principally new features, repre¬ 
senting an interesting topic for conceptual understanding of chaotic dynamics 
and applications. 

In conclusion, let us dwell on open problems. Conhguration space for 
multi-well potentials breaks up naturally into different regions (separate local 
minima); in each of them different dynamical regimes are realized. Thus, 
there is a need for a calculating scheme that will enable us to use the partial 
information about each isolated region with the aim of obtaining a solution 
of the full problem. We believe that the path decomposition expansion of 
Auerbach and Kivelson |in3j will prove to be useful for describing the mixed 
state. This formalism allows the expression of the full time evolution operator 
as a time convolution and surface integrations of products of restricted Green 
functions, each of which involves the sum over paths that are limited to 
different regions of conhguration space. Even for complicated nonseparable 
potentials qualitative behavior is readily inferred and quantitative solutions 
can be obtained from knowledge of the classical dynamics. 

In the late 1990s Zaslavsky and Edelman |104] considered a model of a 
billiard-type system, which consists of two chambers connected through a 
hole. One chamber has a circle-shaped scatterer inside, and the other one 
has a Cassini oval with a concave border. As was shown, the corresponding 
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distribution function does not reach equilibrium even during an anomalous 
large time. We want to emphasize, that the mixed state, at energy a little 
exceeding the saddle one, can serve as a more realistic model for studies of 
anomalous kinetics. 
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